R Code for GCB-19-0810
rm(list = ls()) # Remove all objects from memory
library(ggplot2)
library(gridExtra)
ENSO_colours<-c("blue", "red") # Type of anomaly
my_colours<-c("#feb24c93", "#fb6a4a93", "#9e9ac893", "#a6d96a93") # Thermal regime
theme_set (theme_classic() + theme(panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
axis.line = element_line(colour = "black"),
legend.position="none",
axis.text.x = element_text(angle = 90, vjust = 0.5),
plot.title = element_text(size=12, face="bold"),
#panel.border = element_rect(colour = "black", fill=NA, size=1)
panel.border = element_blank()
))
# Plot ENOS 1971 - 2014
ENOS<-read.table(file="http://www.cpc.ncep.noaa.gov/data/indices/ersst5.nino.mth.81-10.ascii",header=TRUE)
ENOS1982_2014 <-subset(ENOS, YR>1968 & YR<2015) # subset from 1970-2014
ENOS1982_2014$Date<-paste(ENOS1982_2014$YR, ENOS1982_2014$MON, sep = "-" )
ENOS1982_2014$Date<-paste(ENOS1982_2014$Date, "15", sep = "-" )
ENOS1982_2014$Date<-as.Date(ENOS1982_2014$Date, format = "%Y-%m-%d")
ENOS1982_2014$Anomaly <- factor(ifelse(ENOS1982_2014$ANOM.3 >=0, "Positive", "Negative"))
interp <- approx(ENOS1982_2014$Date, ENOS1982_2014$ANOM.3, n=10000)
DataENSO <- data.frame(Date=interp$x, Anomaly3=interp$y)
DataENSO$Anomaly <- factor(ifelse(DataENSO$Anomaly3>=0, "Positive", "Negative"))
ENSO<-ggplot(data=DataENSO, aes(x=Date, y=Anomaly3))+
scale_x_continuous("",expand=c(0, 0))+
scale_y_continuous("SST anomaly (ºC) in Niño.3",expand=c(0, 0))+
geom_area(aes(fill=Anomaly), alpha=0.7)+
scale_fill_manual(values=ENSO_colours) + labs(title="a")+
geom_hline(yintercept=0, linetype="solid", color="black")+
geom_line(position = "identity", size=0.3)+
theme(axis.text.x=element_blank(),
axis.title.x=element_blank(),
axis.ticks.x = element_blank())+
annotate("text", x = 1900, y = 2.0, label = "El Niño \n72-73", size=3)+
annotate("text", x = 5600, y = 2.0, label = "El Niño \n82-83", size=3)+
annotate("text", x = 11100, y = 2.0, label = "El Niño \n97-98", size=3)
ENSO
# ENSOb<-ggplot(data=ENOS1982_2014, aes(x=Date, y=ANOM.3))+
# geom_hline(yintercept=0, linetype="solid", color="black")+
# scale_x_date("", limit=c(as.Date("1969-01-15"), as.Date("2014-12-31")),
# breaks= "12 months",
# expand=c(0, 0))+
# scale_y_continuous("SST anomaly (ºC) in Niño.3",expand=c(0, 0))+
# geom_area(aes(fill=Anomaly))+
# geom_line(position = "identity")+
# scale_fill_manual(values=ENSO_colours) + labs(title="(a)")
# ENSOb
# Select data
cover <- read.csv("Data/Coral_Cover_ETP.csv", header=TRUE) # Select file: Coral_Cover_ETP.csv
# Convert to numeric variable the response variable and convert to nominal variables the factors
cover$Coral_cover <- as.numeric(cover$Coral_cover) # convert to numeric variable
cover$Region <- as.factor(cover$Region) # convert to nominal factors
cover$Location <- as.factor(cover$Location)
cover$Site <- as.factor(cover$Site)
cover$Country <- as.factor(cover$Country)
cover$Thermal_regime <- factor(cover$Thermal_regime,
levels=c("Thermally_Stable", "Upwelling", "Galapagos_Islands","Seasonal"),
labels=c("Thermally stable", "Tropical upwelling", "Equatorial upwelling","Seasonal"))
cover$Lat <- as.numeric(cover$Lat)
cover$Lon <- as.numeric(cover$Lon)
cover$Year <- as.numeric(cover$Year)
#cover$Year <- as.numeric((cover$Year)-1970) # transform years to begin in 0
summary(cover)
## ID1 Reference Province
## Min. : 1.0 Fong, 2017 : 33 Cortezian : 10
## 1st Qu.:144.2 Glynn, 1983 : 31 Galapagos : 61
## Median :285.5 Calderón-Aguilera, 2006: 26 Tropical_East_Pacific:495
## Mean :285.1 Glynn, 1990 : 24
## 3rd Qu.:426.8 Guzman, 2004 : 24
## Max. :568.0 Cortés, 2010 : 23
## (Other) :405
## Ecoregion Region
## Nicoya :238 Gulf_of_Chiriqui :133
## Panama_Bight :117 Colombia : 75
## Cortezian : 71 Mexican_Central : 72
## Mexican_Tropical_Pacific : 50 Costa_Rica_Southern : 56
## Eastern_Galapagos_Islands : 31 Nicaragua_Papagayo_Zone: 49
## Northern_Galapagos_Islands: 28 Panama_Bight : 42
## (Other) : 31 (Other) :139
## Country Platform Thermal_regime
## Panama :151 Continental:470 Thermally stable :250
## Costa_Rica:143 Oceanic : 96 Tropical upwelling :177
## Mexico :123 Equatorial upwelling: 61
## Colombia : 75 Seasonal : 78
## Ecuador : 61
## Nicaragua : 11
## (Other) : 2
## Location Site Year Date
## Gulf_of_Chiriqui:121 Uva_Island : 67 Min. :1970 1/1/09 : 55
## Gorgona_Island : 56 La_Azufrada : 29 1st Qu.:1987 1/1/06 : 32
## Pearl_Islands : 40 Cano_island : 27 Median :1999 1/1/02 : 31
## Bahia_Banderas : 33 Gorgona_Island: 26 Mean :1997 1/1/03 : 31
## Gulf_of_Papagayo: 30 Saboga_Island : 24 3rd Qu.:2006 1/1/97 : 26
## Los_Cabos : 28 Cabo_Pulmo : 17 Max. :2014 1/1/98 : 22
## (Other) :258 (Other) :376 (Other):369
## Period_ENOS Coral_cover Algal_cover Lat
## 1971-1972: 24 Min. : 0.000 Min. : 0.00 Min. :-1.352
## 1973-1982: 60 1st Qu.: 8.568 1st Qu.:12.61 1st Qu.: 5.551
## 1983-1997:169 Median :25.000 Median :35.84 Median : 8.578
## 1998-2016:313 Mean :29.185 Mean :37.57 Mean : 9.499
## 3rd Qu.:45.725 3rd Qu.:54.81 3rd Qu.:10.791
## Max. :95.720 Max. :98.00 Max. :25.796
## NA's :390
## Lon
## Min. :-111.23
## 1st Qu.: -91.82
## Median : -83.87
## Mean : -87.92
## 3rd Qu.: -81.75
## Max. : -77.34
##
# Aggregate by site
aggr.site <- aggregate(Coral_cover ~ Year+Region+Location+Site+Period_ENOS+Thermal_regime,FUN=mean,data=cover)
# Aggregate by location
aggr.location <- aggregate(Coral_cover~Year+Region+Location+Period_ENOS+Thermal_regime,FUN=mean,data=aggr.site)
# Aggregate by region
aggr.region <- aggregate(Coral_cover~Year+Region+Period_ENOS+Thermal_regime,FUN=mean,data=aggr.location)
# Create groups for each of three ENSO time-intervals
ENSO.73_82 <- subset(aggr.region, Period_ENOS == "1973-1982")
ENSO.83_97 <- subset(aggr.region, Period_ENOS == "1983-1997")
ENSO.98_14 <- subset(aggr.region, Period_ENOS == "1998-2016")
# Exclude regions with cero and only one observation of coral cover per time-interval
ENSO.83_97<- ENSO.83_97[!ENSO.83_97$Region == "Costa_Rica_Central", ]
ENSO.83_97<- ENSO.83_97[!ENSO.83_97$Region == "El_Salvador", ]
ENSO.83_97<- ENSO.83_97[!ENSO.83_97$Region == "Western_Galapagos_Islands", ]
ENSO.83_97<- ENSO.83_97[!ENSO.83_97$Region == "Clipperton", ]
ENSO.83_97<- ENSO.83_97[!ENSO.83_97$Region == "Mexican_Southern", ]
ENSO.83_97<- ENSO.83_97[!ENSO.83_97$Region == "Northern_Galapagos_Islands", ]
ENSO.98_14<- ENSO.98_14[!ENSO.98_14$Region == "Clipperton", ]
ENSO.98_14<- ENSO.98_14[!ENSO.98_14$Region == "El_Salvador", ]
ENSO.98_14<- ENSO.98_14[!ENSO.98_14$Region == "Western_Galapagos_Islands", ]
ENSO.98_14<- ENSO.98_14[!ENSO.98_14$Region == "Revillagigedo_Islands", ]
ENSO.83_14 <- rbind(ENSO.83_97,ENSO.98_14) # Second and third time intervals
ENSO.73_14 <- rbind(ENSO.73_82,ENSO.83_14) # First, second and third time intervals
ENSO.73_97 <- rbind(ENSO.73_82,ENSO.83_97) # First and second and time intervals
library(nlme)
# All years 1970_2014
model1 <- lme(
Coral_cover ~ -1 + Thermal_regime * Year, random = ~1|Region, data=aggr.region)
summary(model1)
## Linear mixed-effects model fit by REML
## Data: aggr.region
## AIC BIC logLik
## 1727.902 1760.58 -853.9508
##
## Random effects:
## Formula: ~1 | Region
## (Intercept) Residual
## StdDev: 12.5915 16.52254
##
## Fixed effects: Coral_cover ~ -1 + Thermal_regime * Year
## Value Std.Error DF t-value
## Thermal_regimeThermally stable -428.2896 377.6124 181 -1.1342042
## Thermal_regimeTropical upwelling -323.5241 371.8720 14 -0.8699878
## Thermal_regimeEquatorial upwelling 254.5442 573.7481 14 0.4436514
## Thermal_regimeSeasonal 1036.7134 843.3488 181 1.2292819
## Year 0.2264 0.1889 181 1.1987453
## Thermal_regimeTropical upwelling:Year -0.0478 0.2650 181 -0.1804381
## Thermal_regimeEquatorial upwelling:Year -0.3448 0.3436 181 -1.0033294
## Thermal_regimeSeasonal:Year -0.7293 0.4592 181 -1.5880369
## p-value
## Thermal_regimeThermally stable 0.2582
## Thermal_regimeTropical upwelling 0.3990
## Thermal_regimeEquatorial upwelling 0.6641
## Thermal_regimeSeasonal 0.2206
## Year 0.2322
## Thermal_regimeTropical upwelling:Year 0.8570
## Thermal_regimeEquatorial upwelling:Year 0.3170
## Thermal_regimeSeasonal:Year 0.1140
## Correlation:
## Thr_Ts Thr_Tu Thr_Eu Thrm_S Year
## Thermal_regimeTropical upwelling 0.000
## Thermal_regimeEquatorial upwelling 0.000 0.000
## Thermal_regimeSeasonal 0.014 0.000 0.000
## Year -1.000 0.000 0.000 -0.015
## Thermal_regimeTropical upwelling:Year 0.713 -0.701 0.000 0.010 -0.713
## Thermal_regimeEquatorial upwelling:Year 0.550 0.000 -0.835 0.008 -0.550
## Thermal_regimeSeasonal:Year 0.398 0.000 0.000 -0.912 -0.397
## T_Tu:Y T_Eu:Y
## Thermal_regimeTropical upwelling
## Thermal_regimeEquatorial upwelling
## Thermal_regimeSeasonal
## Year
## Thermal_regimeTropical upwelling:Year
## Thermal_regimeEquatorial upwelling:Year 0.392
## Thermal_regimeSeasonal:Year 0.283 0.218
##
## Standardized Within-Group Residuals:
## Min Q1 Med Q3 Max
## -2.323212218 -0.685668507 -0.007907488 0.578042210 3.071703318
##
## Number of Observations: 202
## Number of Groups: 16
anova(model1)
## numDF denDF F-value p-value
## Thermal_regime 4 179 14.698982 <.0001
## Year 1 179 0.717103 0.3982
## Thermal_regime:Year 3 179 1.093597 0.3532
plot(ranef(model1)) # shows symmetrical scatter effects around zero
plot(model1) # plot residuals vs fitted
resnorm1 <- resid(model1)
hist(resnorm1, xlab = "Residuals", main = "") # residuals are normally distributed?
coef.m1 <- as.data.frame(coef(summary(model1))) # Coefficients of the model
plot(model1, resid(., type = "p") ~ fitted(.) | Region, abline = 0)
plot(model1, Coral_cover ~ fitted(.) | Region, abline = c(0,1))
plot(model1, Coral_cover ~ fitted(.) | Thermal_regime, abline = c(0,1))
plot(model1, Coral_cover ~ fitted(.), abline = c(0,1))
#Using raw data (Only thermal regime is significant in model)
Model1_plot <- ggplot(aggr.region, aes(x=Year, y=Coral_cover)) +
geom_boxplot(aes(group=Year), outlier.shape = NA) +
stat_boxplot(aes(group=Year), geom = 'errorbar')+
geom_smooth(method = lm, se=FALSE, linetype = "dashed", colour="red")+
geom_smooth(span = 0.2, se=T, colour="darkgray")+
#geom_point(aes(colour=Thermal_regime), alpha=0.5) + ylim(0,79) +
scale_y_continuous("Live coral cover (%)", limits = c(-2, 80), expand = c(0,0))+
scale_x_continuous("", limits = c(1969, 2015),
breaks = seq(1970, 2014, by=2), expand = c(0,0))+
scale_color_manual(values=my_colours) + labs(title="b")+
#annotate("text", x = 1983, y = 70, label = "*", size=8)+
annotate("rect", xmin = 1982, xmax = 1983, ymin = 0, ymax = 80, alpha = .2, fill="red")+
annotate("rect", xmin = 1997, xmax = 1998, ymin = 0, ymax = 80, alpha = .2, fill="red")+
annotate("rect", xmin = 1972, xmax = 1973, ymin = 0, ymax = 80, alpha = .2, fill="red")
Model1_plot
# Time-interval 1 - ENSO.73_82
model2 <- lme(Coral_cover ~ -1 + Thermal_regime * Year, random = ~1|Region, data=ENSO.73_82)
summary(model2)
## Linear mixed-effects model fit by REML
## Data: ENSO.73_82
## AIC BIC logLik
## 164.3022 170.4829 -74.15108
##
## Random effects:
## Formula: ~1 | Region
## (Intercept) Residual
## StdDev: 20.41219 11.23257
##
## Fixed effects: Coral_cover ~ -1 + Thermal_regime * Year
## Value Std.Error DF t-value
## Thermal_regimeThermally stable -570.074 2795.777 4 -0.2039053
## Thermal_regimeTropical upwelling -8303.810 3542.253 4 -2.3442170
## Thermal_regimeEquatorial upwelling -3307.477 3419.771 4 -0.9671634
## Year 0.301 1.412 13 0.2131857
## Thermal_regimeTropical upwelling:Year 3.912 2.281 13 1.7150859
## Thermal_regimeEquatorial upwelling:Year 1.381 2.232 13 0.6189979
## p-value
## Thermal_regimeThermally stable 0.8484
## Thermal_regimeTropical upwelling 0.0790
## Thermal_regimeEquatorial upwelling 0.3882
## Year 0.8345
## Thermal_regimeTropical upwelling:Year 0.1101
## Thermal_regimeEquatorial upwelling:Year 0.5466
## Correlation:
## Thr_Ts Thr_Tu Thr_Eu Year T_Tu:Y
## Thermal_regimeTropical upwelling 0.000
## Thermal_regimeEquatorial upwelling 0.000 0.000
## Year -1.000 0.000 0.000
## Thermal_regimeTropical upwelling:Year 0.619 -0.785 0.000 -0.619
## Thermal_regimeEquatorial upwelling:Year 0.633 0.000 -0.774 -0.633 0.392
##
## Standardized Within-Group Residuals:
## Min Q1 Med Q3 Max
## -1.4650132 -0.4802879 -0.2103293 0.5465022 1.8906401
##
## Number of Observations: 22
## Number of Groups: 7
anova(model2)
## numDF denDF F-value p-value
## Thermal_regime 3 4 3.434033 0.1322
## Year 1 13 3.583180 0.0808
## Thermal_regime:Year 2 13 1.472433 0.2652
plot(ranef(model2))
plot(model2)
resnorm2 <- resid(model2)
hist(resnorm2, xlab = "Residuals", main = "")
coef.m2 <- as.data.frame(coef(summary(model2)))
plot(model2, resid(., type = "p") ~ fitted(.) | Region, abline = 0)
plot(model2, Coral_cover ~ fitted(.) | Region, abline = c(0,1))
plot(model2, Coral_cover ~ fitted(.) | Thermal_regime, abline = c(0,1))
Note that this first time interval has low number of regions at each year, even the GLMM model is adjusted, please warning with interpretation/generalization of results
#Using raw data (the model is not significant)
Model2_plot <- ggplot(ENSO.73_82, aes(x=Year, y=Coral_cover, colour=Thermal_regime)) +
geom_smooth(method = lm, se=FALSE, linetype = "dashed", alpha=1, size=2)+
geom_point(aes(fill=factor(Thermal_regime)),
shape = 21, colour = "black", size = 2, stroke = 0.3, alpha=0.5) +
scale_y_continuous("Live coral cover (%)", limits = c(0, 80))+
scale_x_continuous("", limits = c(1973.7, 1982.3),
breaks = seq(1974, 1982, by=2), expand = c(0.02,0.02))+
scale_fill_manual(values=my_colours) +
scale_colour_manual(values=my_colours)+
labs(title="c")
Model2_plot
# Time-interval 3 - ENSO.83_97
model3 <- lme(Coral_cover ~ -1 + Thermal_regime * Year, random = ~1|Region, data=ENSO.83_97)
summary(model3)
## Linear mixed-effects model fit by REML
## Data: ENSO.83_97
## AIC BIC logLik
## 476.7179 497.1484 -228.3589
##
## Random effects:
## Formula: ~1 | Region
## (Intercept) Residual
## StdDev: 4.741945 9.686494
##
## Fixed effects: Coral_cover ~ -1 + Thermal_regime * Year
## Value Std.Error DF t-value
## Thermal_regimeThermally stable -2536.003 926.2417 50 -2.737950
## Thermal_regimeTropical upwelling -7100.953 882.5060 8 -8.046351
## Thermal_regimeEquatorial upwelling 111.787 2814.6054 8 0.039717
## Thermal_regimeSeasonal 6083.292 1940.6223 50 3.134712
## Year 1.281 0.4654 50 2.753386
## Thermal_regimeTropical upwelling:Year 2.306 0.6426 50 3.589323
## Thermal_regimeEquatorial upwelling:Year -1.337 1.4905 50 -0.897275
## Thermal_regimeSeasonal:Year -4.319 1.0734 50 -4.023670
## p-value
## Thermal_regimeThermally stable 0.0085
## Thermal_regimeTropical upwelling 0.0000
## Thermal_regimeEquatorial upwelling 0.9693
## Thermal_regimeSeasonal 0.0029
## Year 0.0082
## Thermal_regimeTropical upwelling:Year 0.0008
## Thermal_regimeEquatorial upwelling:Year 0.3739
## Thermal_regimeSeasonal:Year 0.0002
## Correlation:
## Thr_Ts Thr_Tu Thr_Eu Thrm_S Year
## Thermal_regimeTropical upwelling 0.000
## Thermal_regimeEquatorial upwelling 0.000 0.000
## Thermal_regimeSeasonal 0.014 0.000 0.000
## Year -1.000 0.000 0.000 -0.014
## Thermal_regimeTropical upwelling:Year 0.724 -0.689 0.000 0.010 -0.724
## Thermal_regimeEquatorial upwelling:Year 0.312 0.000 -0.950 0.005 -0.312
## Thermal_regimeSeasonal:Year 0.420 0.000 0.000 -0.901 -0.420
## T_Tu:Y T_Eu:Y
## Thermal_regimeTropical upwelling
## Thermal_regimeEquatorial upwelling
## Thermal_regimeSeasonal
## Year
## Thermal_regimeTropical upwelling:Year
## Thermal_regimeEquatorial upwelling:Year 0.226
## Thermal_regimeSeasonal:Year 0.304 0.131
##
## Standardized Within-Group Residuals:
## Min Q1 Med Q3 Max
## -2.99663110 -0.53542294 -0.06265663 0.58765352 2.99408354
##
## Number of Observations: 65
## Number of Groups: 10
anova(model3)
## numDF denDF F-value p-value
## Thermal_regime 4 48 62.03136 <.0001
## Year 1 48 39.45390 <.0001
## Thermal_regime:Year 3 48 14.56275 <.0001
plot(ranef(model3))
plot(model3)
resnorm3 <- resid(model3)
hist(resnorm3, xlab = "Residuals", main = "")
coef.m3 <- as.data.frame(coef(summary(model3)))
plot(model3, resid(., type = "p") ~ fitted(.) | Region, abline = 0)
plot(model3, Coral_cover ~ fitted(.) | Region, abline = c(0,1))
plot(model3, Coral_cover ~ fitted(.) | Thermal_regime, abline = c(0,1))
# Create a new data frame for independent variables
NewData_3 <- expand.grid(Thermal_regime=unique(ENSO.83_97$Thermal_regime),
#Region=unique(aggr.region$Region),
Year=seq((min(ENSO.83_97$Year)), (max(ENSO.83_97$Year))))
pred_3 <- predict(model3, newdata=NewData_3, level=0, asList = FALSE)
# Using model data
Model3b_plot <- ggplot(ENSO.83_97, aes(x=Year, y=Coral_cover, colour=Thermal_regime)) +
geom_point(aes(fill=factor(Thermal_regime)),
shape = 21, colour = "black", size = 2, stroke = 0.3, alpha=0.5) +
scale_y_continuous(" ", limits = c(0, 80))+
scale_x_continuous("", limits = c(1982.7, 1997.3),
breaks = seq(1983, 1997, by=2), expand = c(0.02,0.02))+
geom_line(data=NewData_3,
aes(y=predict(model3, level=0, newdata=NewData_3)), size=2) +
theme(axis.text.y=element_blank(),
axis.title.y=element_blank(),
legend.position = c(0.58, -0.4),
legend.title=element_blank(),
legend.text=element_text(size=6),
legend.box ="horizontal")+
guides(col = guide_legend(nrow = 2))+
scale_fill_manual(values=my_colours) +
scale_colour_manual(values=my_colours)+
labs(title="d")
Model3b_plot
# Time-interval 4 - ENSO.98_14
model4 <- lme(Coral_cover ~ -1 + Thermal_regime * Year, random = ~1|Region, data=ENSO.98_14)
summary(model4)
## Linear mixed-effects model fit by REML
## Data: ENSO.98_14
## AIC BIC logLik
## 866.8707 892.5142 -423.4354
##
## Random effects:
## Formula: ~1 | Region
## (Intercept) Residual
## StdDev: 11.83418 15.29875
##
## Fixed effects: Coral_cover ~ -1 + Thermal_regime * Year
## Value Std.Error DF t-value
## Thermal_regimeThermally stable -1881.9759 1396.7590 87 -1.3473877
## Thermal_regimeTropical upwelling 2994.9557 1056.7535 10 2.8341102
## Thermal_regimeEquatorial upwelling -1881.0600 2875.6287 10 -0.6541387
## Thermal_regimeSeasonal 247.4085 1577.1131 87 0.1568743
## Year 0.9516 0.6965 87 1.3661224
## Thermal_regimeTropical upwelling:Year -2.4271 0.8732 87 -2.7795278
## Thermal_regimeEquatorial upwelling:Year -0.0018 1.5937 87 -0.0011048
## Thermal_regimeSeasonal:Year -1.0581 1.0458 87 -1.0117420
## p-value
## Thermal_regimeThermally stable 0.1814
## Thermal_regimeTropical upwelling 0.0177
## Thermal_regimeEquatorial upwelling 0.5278
## Thermal_regimeSeasonal 0.8757
## Year 0.1754
## Thermal_regimeTropical upwelling:Year 0.0067
## Thermal_regimeEquatorial upwelling:Year 0.9991
## Thermal_regimeSeasonal:Year 0.3145
## Correlation:
## Thr_Ts Thr_Tu Thr_Eu Thrm_S Year
## Thermal_regimeTropical upwelling 0.000
## Thermal_regimeEquatorial upwelling 0.000 0.000
## Thermal_regimeSeasonal 0.008 0.000 0.000
## Year -1.000 0.000 0.000 -0.008
## Thermal_regimeTropical upwelling:Year 0.798 -0.603 0.000 0.006 -0.798
## Thermal_regimeEquatorial upwelling:Year 0.437 0.000 -0.899 0.003 -0.437
## Thermal_regimeSeasonal:Year 0.660 0.000 0.000 -0.746 -0.660
## T_Tu:Y T_Eu:Y
## Thermal_regimeTropical upwelling
## Thermal_regimeEquatorial upwelling
## Thermal_regimeSeasonal
## Year
## Thermal_regimeTropical upwelling:Year
## Thermal_regimeEquatorial upwelling:Year 0.349
## Thermal_regimeSeasonal:Year 0.526 0.288
##
## Standardized Within-Group Residuals:
## Min Q1 Med Q3 Max
## -1.98327607 -0.50110676 -0.02557425 0.30976474 2.71954447
##
## Number of Observations: 104
## Number of Groups: 12
anova(model4)
## numDF denDF F-value p-value
## Thermal_regime 4 85 15.541089 <.0001
## Year 1 85 1.227329 0.2711
## Thermal_regime:Year 3 85 2.983540 0.0358
plot(ranef(model4))
plot(model4)
resnorm4 <- resid(model4)
hist(resnorm4, xlab = "Residuals", main = "")
coef.m4 <- as.data.frame(coef(summary(model4)))
plot(model4, resid(., type = "p") ~ fitted(.) | Region, abline = 0)
plot(model4, Coral_cover ~ fitted(.) | Region, abline = c(0,1))
plot(model4, Coral_cover ~ fitted(.) | Thermal_regime, abline = c(0,1))
# Competing model with autocorrelation structure
model4b <- lme(Coral_cover ~ -1 + Thermal_regime * Year, random = ~1|Region,
data=ENSO.98_14, correlation=corCompSymm(form=~Year|Region))
anova(model4, model4b) # ~ same results
## Model df AIC BIC logLik Test L.Ratio p-value
## model4 1 10 866.8707 892.5142 -423.4354
## model4b 2 11 868.8707 897.0786 -423.4354 1 vs 2 4.547474e-13 1
# Create a new data frame for independent variables
NewData_4 <- expand.grid(Thermal_regime=unique(ENSO.98_14$Thermal_regime),
#Region=unique(aggr.region$Region),
Year=seq((min(ENSO.98_14$Year)), (max(ENSO.98_14$Year))))
pred_4 <- predict(model4, newdata=NewData_4, level=0)
#summary(pred)
#length(pred)
# Using model data
Model4b_plot <- ggplot(ENSO.98_14, aes(x=Year, y=Coral_cover, colour=Thermal_regime)) +
geom_point(aes(fill=factor(Thermal_regime)),
shape = 21, colour = "black", size = 2, stroke = 0.3, alpha=0.5) +
scale_y_continuous(" ", limits = c(0, 80))+
scale_x_continuous("", limits = c(1997.7, 2014.3),
breaks = seq(1998, 2014, by=2), expand = c(0.02,0.02))+
geom_line(data=NewData_4, aes(y=predict(model4, level=0, newdata=NewData_4)), size=2)+
theme(axis.text.y=element_blank(), axis.title.y=element_blank())+
scale_fill_manual(values=my_colours) +
scale_colour_manual(values=my_colours)+
labs(title="e")
Model4b_plot
# Time-interval 5 - ENSO.83_14
model5 <- lme(Coral_cover ~ -1 + Thermal_regime * Year, random = ~1|Region, data=ENSO.83_14)
summary(model5)
## Linear mixed-effects model fit by REML
## Data: ENSO.83_14
## AIC BIC logLik
## 1423.169 1453.983 -701.5844
##
## Random effects:
## Formula: ~1 | Region
## (Intercept) Residual
## StdDev: 12.53351 15.47206
##
## Fixed effects: Coral_cover ~ -1 + Thermal_regime * Year
## Value Std.Error DF t-value
## Thermal_regimeThermally stable -1736.8288 514.4829 151 -3.375873
## Thermal_regimeTropical upwelling 36.1048 507.2122 11 0.071183
## Thermal_regimeEquatorial upwelling -1334.0770 907.7813 11 -1.469602
## Thermal_regimeSeasonal 1031.1175 792.1371 151 1.301691
## Year 0.8796 0.2572 151 3.420336
## Thermal_regimeTropical upwelling:Year -0.8804 0.3609 151 -2.439274
## Thermal_regimeEquatorial upwelling:Year -0.2024 0.5212 151 -0.388367
## Thermal_regimeSeasonal:Year -1.3793 0.4685 151 -2.943936
## p-value
## Thermal_regimeThermally stable 0.0009
## Thermal_regimeTropical upwelling 0.9445
## Thermal_regimeEquatorial upwelling 0.1697
## Thermal_regimeSeasonal 0.1950
## Year 0.0008
## Thermal_regimeTropical upwelling:Year 0.0159
## Thermal_regimeEquatorial upwelling:Year 0.6983
## Thermal_regimeSeasonal:Year 0.0038
## Correlation:
## Thr_Ts Thr_Tu Thr_Eu Thrm_S Year
## Thermal_regimeTropical upwelling 0.000
## Thermal_regimeEquatorial upwelling 0.000 0.000
## Thermal_regimeSeasonal 0.015 0.000 0.000
## Year -1.000 0.000 0.000 -0.016
## Thermal_regimeTropical upwelling:Year 0.712 -0.702 0.000 0.011 -0.712
## Thermal_regimeEquatorial upwelling:Year 0.493 0.000 -0.870 0.008 -0.493
## Thermal_regimeSeasonal:Year 0.536 0.000 0.000 -0.836 -0.535
## T_Tu:Y T_Eu:Y
## Thermal_regimeTropical upwelling
## Thermal_regimeEquatorial upwelling
## Thermal_regimeSeasonal
## Year
## Thermal_regimeTropical upwelling:Year
## Thermal_regimeEquatorial upwelling:Year 0.352
## Thermal_regimeSeasonal:Year 0.381 0.264
##
## Standardized Within-Group Residuals:
## Min Q1 Med Q3 Max
## -2.52794002 -0.47964749 -0.06607999 0.50350339 3.04770181
##
## Number of Observations: 169
## Number of Groups: 13
anova(model5)
## numDF denDF F-value p-value
## Thermal_regime 4 149 14.198648 <.0001
## Year 1 149 4.270123 0.0405
## Thermal_regime:Year 3 149 3.799303 0.0116
plot(ranef(model5))
plot(model5)
resnorm5 <- resid(model5)
hist(resnorm5, xlab = "Residuals", main = "")
coef.m5 <- as.data.frame(coef(summary(model5)))
plot(model5, resid(., type = "p") ~ fitted(.) | Region, abline = 0)
plot(model5, Coral_cover ~ fitted(.) | Region, abline = c(0,1))
plot(model5, Coral_cover ~ fitted(.) | Thermal_regime, abline = c(0,1))
# Create a new data frame for independent variables
NewData_5 <- expand.grid(Thermal_regime=unique(ENSO.83_14$Thermal_regime),
#Region=unique(aggr.region$Region),
Year=seq((min(ENSO.83_14$Year)), (max(ENSO.83_14$Year))))
pred_5 <- predict(model5, newdata=NewData_5, level=0) # level zero=population predictions
# Using model data
Model5b_plot <- ggplot(ENSO.83_14, aes(x=Year, y=Coral_cover, colour=Thermal_regime)) +
geom_point(aes(fill=factor(Thermal_regime)),
shape = 21, colour = "black", size = 2, stroke = 0.3, alpha=0.5) +
scale_y_continuous("Live coral cover (%)", limits = c(0, 80))+
scale_x_continuous("", limits = c(1982.7, 2014.3),
breaks = seq(1983, 2014, by=2), expand = c(0.02,0.02))+
geom_line(data=NewData_5, aes(y=predict(model5, level=0, newdata=NewData_5)), size=2)+
theme(axis.text.y=element_blank(),
axis.title.y=element_blank(),
legend.position = c(0.58, -0.3),
legend.title=element_blank(),
legend.text=element_text(size=6),
legend.direction="horizontal")+
scale_fill_manual(values=my_colours) +
scale_colour_manual(values=my_colours)+
labs(title="Model 5")
Model5b_plot
library(multcomp)
library(lme4)
ENSO.73_97$Year_F <- as.factor(ENSO.73_97$Year)
comp_years_ENSO.73_97 <- aov(Coral_cover ~ Year_F + Thermal_regime, data = ENSO.73_97)
m.comp_ENSO.73_97 <- glht(comp_years_ENSO.73_97, linfct = mcp(Year_F = "Tukey"))
summary(m.comp_ENSO.73_97, test = univariate())
##
## Simultaneous Tests for General Linear Hypotheses
##
## Multiple Comparisons of Means: Tukey Contrasts
##
##
## Fit: aov(formula = Coral_cover ~ Year_F + Thermal_regime, data = ENSO.73_97)
##
## Linear Hypotheses:
## Estimate Std. Error t value Pr(>|t|)
## 1975 - 1974 == 0 -3.44519 12.08474 -0.285 0.776577
## 1976 - 1974 == 0 11.67379 13.58569 0.859 0.393669
## 1977 - 1974 == 0 2.20768 17.29237 0.128 0.898846
## 1978 - 1974 == 0 -21.04773 13.64079 -1.543 0.128178
## 1979 - 1974 == 0 12.32727 13.64079 0.904 0.369826
## 1980 - 1974 == 0 2.61704 12.23990 0.214 0.831431
## 1981 - 1974 == 0 10.89268 17.29237 0.630 0.531184
## 1982 - 1974 == 0 9.29041 10.82386 0.858 0.394186
## 1983 - 1974 == 0 -26.56453 12.25535 -2.168 0.034237 *
## 1984 - 1974 == 0 -17.28488 11.32342 -1.526 0.132236
## 1985 - 1974 == 0 -6.66044 13.76573 -0.484 0.630289
## 1986 - 1974 == 0 -16.73332 11.32342 -1.478 0.144790
## 1987 - 1974 == 0 -6.97546 10.29300 -0.678 0.500618
## 1988 - 1974 == 0 -5.21764 11.53808 -0.452 0.652776
## 1989 - 1974 == 0 -9.53355 10.53984 -0.905 0.369396
## 1990 - 1974 == 0 -10.69232 17.29237 -0.618 0.538741
## 1991 - 1974 == 0 4.34819 12.74709 0.341 0.734231
## 1992 - 1974 == 0 0.52130 10.65659 0.049 0.961150
## 1993 - 1974 == 0 -6.10824 11.04278 -0.553 0.582255
## 1994 - 1974 == 0 3.00209 11.00761 0.273 0.786014
## 1995 - 1974 == 0 25.24915 11.51467 2.193 0.032279 *
## 1996 - 1974 == 0 8.32457 11.51467 0.723 0.472565
## 1997 - 1974 == 0 0.07463 10.75923 0.007 0.994489
## 1998 - 1974 == 0 26.29686 17.32517 1.518 0.134395
## 1976 - 1975 == 0 15.11898 13.58569 1.113 0.270282
## 1977 - 1975 == 0 5.65287 17.29237 0.327 0.744901
## 1978 - 1975 == 0 -17.60254 13.64079 -1.290 0.201934
## 1979 - 1975 == 0 15.77246 13.64079 1.156 0.252231
## 1980 - 1975 == 0 6.06223 12.23990 0.495 0.622239
## 1981 - 1975 == 0 14.33787 17.29237 0.829 0.410365
## 1982 - 1975 == 0 12.73559 10.82386 1.177 0.244072
## 1983 - 1975 == 0 -23.11934 12.25535 -1.886 0.064156 .
## 1984 - 1975 == 0 -13.83970 11.32342 -1.222 0.226485
## 1985 - 1975 == 0 -3.21526 13.76573 -0.234 0.816128
## 1986 - 1975 == 0 -13.28814 11.32342 -1.174 0.245307
## 1987 - 1975 == 0 -3.53027 10.29300 -0.343 0.732834
## 1988 - 1975 == 0 -1.77246 11.53808 -0.154 0.878435
## 1989 - 1975 == 0 -6.08837 10.53984 -0.578 0.565697
## 1990 - 1975 == 0 -7.24713 17.29237 -0.419 0.676669
## 1991 - 1975 == 0 7.79338 12.74709 0.611 0.543293
## 1992 - 1975 == 0 3.96649 10.65659 0.372 0.711070
## 1993 - 1975 == 0 -2.66305 11.04278 -0.241 0.810269
## 1994 - 1975 == 0 6.44728 11.00761 0.586 0.560304
## 1995 - 1975 == 0 28.69434 11.51467 2.492 0.015533 *
## 1996 - 1975 == 0 11.76975 11.51467 1.022 0.310882
## 1997 - 1975 == 0 3.51982 10.75923 0.327 0.744717
## 1998 - 1975 == 0 29.74205 17.32517 1.717 0.091281 .
## 1977 - 1976 == 0 -9.46611 18.35985 -0.516 0.608069
## 1978 - 1976 == 0 -32.72152 15.11308 -2.165 0.034435 *
## 1979 - 1976 == 0 0.65348 15.11308 0.043 0.965657
## 1980 - 1976 == 0 -9.05675 13.81029 -0.656 0.514504
## 1981 - 1976 == 0 -0.78111 18.35985 -0.043 0.966209
## 1982 - 1976 == 0 -2.38339 12.41250 -0.192 0.848390
## 1983 - 1976 == 0 -38.23832 13.92671 -2.746 0.007992 **
## 1984 - 1976 == 0 -28.95867 12.90860 -2.243 0.028641 *
## 1985 - 1976 == 0 -18.33423 15.08488 -1.215 0.229053
## 1986 - 1976 == 0 -28.40711 12.90860 -2.201 0.031688 *
## 1987 - 1976 == 0 -18.64925 12.08019 -1.544 0.127987
## 1988 - 1976 == 0 -16.89143 13.27936 -1.272 0.208362
## 1989 - 1976 == 0 -21.20734 12.30424 -1.724 0.090021 .
## 1990 - 1976 == 0 -22.36611 18.35985 -1.218 0.227994
## 1991 - 1976 == 0 -7.32560 14.34376 -0.511 0.611455
## 1992 - 1976 == 0 -11.15249 12.45224 -0.896 0.374095
## 1993 - 1976 == 0 -17.78203 12.71582 -1.398 0.167222
## 1994 - 1976 == 0 -8.67170 12.78086 -0.678 0.500113
## 1995 - 1976 == 0 13.57536 13.30638 1.020 0.311792
## 1996 - 1976 == 0 -3.34922 13.30638 -0.252 0.802147
## 1997 - 1976 == 0 -11.59916 12.56356 -0.923 0.359644
## 1998 - 1976 == 0 14.62307 18.62181 0.785 0.435440
## 1978 - 1977 == 0 -23.25541 18.23689 -1.275 0.207242
## 1979 - 1977 == 0 10.11959 18.23689 0.555 0.581063
## 1980 - 1977 == 0 0.40936 17.14223 0.024 0.981029
## 1981 - 1977 == 0 8.68500 20.93138 0.415 0.679700
## 1982 - 1977 == 0 7.08272 16.43022 0.431 0.667982
## 1983 - 1977 == 0 -28.77221 17.29679 -1.663 0.101527
## 1984 - 1977 == 0 -19.49257 16.66532 -1.170 0.246846
## 1985 - 1977 == 0 -8.86813 18.12710 -0.489 0.626500
## 1986 - 1977 == 0 -18.94101 16.66532 -1.137 0.260321
## 1987 - 1977 == 0 -9.18314 15.95407 -0.576 0.567075
## 1988 - 1977 == 0 -7.42533 16.79651 -0.442 0.660051
## 1989 - 1977 == 0 -11.74124 16.16342 -0.726 0.470461
## 1990 - 1977 == 0 -12.90000 20.93138 -0.616 0.540067
## 1991 - 1977 == 0 2.14051 17.69601 0.121 0.904134
## 1992 - 1977 == 0 -1.68638 16.10151 -0.105 0.916942
## 1993 - 1977 == 0 -8.31592 16.53313 -0.503 0.616847
## 1994 - 1977 == 0 0.79441 16.37629 0.049 0.961474
## 1995 - 1977 == 0 23.04147 16.81692 1.370 0.175836
## 1996 - 1977 == 0 6.11688 16.81692 0.364 0.717357
## 1997 - 1977 == 0 -2.13305 16.22072 -0.132 0.895826
## 1998 - 1977 == 0 24.08918 21.30941 1.130 0.262862
## 1979 - 1978 == 0 33.37500 14.80072 2.255 0.027861 *
## 1980 - 1978 == 0 23.66477 13.52755 1.749 0.085426 .
## 1981 - 1978 == 0 31.94041 18.23689 1.751 0.085071 .
## 1982 - 1978 == 0 30.33813 12.59792 2.408 0.019176 *
## 1983 - 1978 == 0 -5.51680 13.52755 -0.408 0.684882
## 1984 - 1978 == 0 3.76284 12.90860 0.291 0.771692
## 1985 - 1978 == 0 14.38729 14.93497 0.963 0.339314
## 1986 - 1978 == 0 4.31440 12.90860 0.334 0.739393
## 1987 - 1978 == 0 14.07227 11.92988 1.180 0.242901
## 1988 - 1978 == 0 15.83008 12.90266 1.227 0.224739
## 1989 - 1978 == 0 11.51417 12.16281 0.947 0.347668
## 1990 - 1978 == 0 10.35541 18.23689 0.568 0.572307
## 1991 - 1978 == 0 25.39592 14.02252 1.811 0.075220 .
## 1992 - 1978 == 0 21.56903 12.12366 1.779 0.080377 .
## 1993 - 1978 == 0 14.93949 12.64762 1.181 0.242259
## 1994 - 1978 == 0 24.04982 12.43222 1.934 0.057855 .
## 1995 - 1978 == 0 46.29688 12.85667 3.601 0.000651 ***
## 1996 - 1978 == 0 29.37229 12.85667 2.285 0.025949 *
## 1997 - 1978 == 0 21.12236 12.22384 1.728 0.089226 .
## 1998 - 1978 == 0 47.34459 18.23689 2.596 0.011882 *
## 1980 - 1979 == 0 -9.71023 13.52755 -0.718 0.475706
## 1981 - 1979 == 0 -1.43459 18.23689 -0.079 0.937566
## 1982 - 1979 == 0 -3.03687 12.59792 -0.241 0.810344
## 1983 - 1979 == 0 -38.89180 13.52755 -2.875 0.005611 **
## 1984 - 1979 == 0 -29.61216 12.90860 -2.294 0.025369 *
## 1985 - 1979 == 0 -18.98771 14.93497 -1.271 0.208590
## 1986 - 1979 == 0 -29.06060 12.90860 -2.251 0.028108 *
## 1987 - 1979 == 0 -19.30273 11.92988 -1.618 0.110993
## 1988 - 1979 == 0 -17.54492 12.90266 -1.360 0.179070
## 1989 - 1979 == 0 -21.86083 12.16281 -1.797 0.077398 .
## 1990 - 1979 == 0 -23.01959 18.23689 -1.262 0.211824
## 1991 - 1979 == 0 -7.97908 14.02252 -0.569 0.571503
## 1992 - 1979 == 0 -11.80597 12.12366 -0.974 0.334134
## 1993 - 1979 == 0 -18.43551 12.64762 -1.458 0.150246
## 1994 - 1979 == 0 -9.32518 12.43222 -0.750 0.456186
## 1995 - 1979 == 0 12.92188 12.85667 1.005 0.318967
## 1996 - 1979 == 0 -4.00271 12.85667 -0.311 0.756645
## 1997 - 1979 == 0 -12.25264 12.22384 -1.002 0.320265
## 1998 - 1979 == 0 13.96959 18.23689 0.766 0.446726
## 1981 - 1980 == 0 8.27564 17.14223 0.483 0.631050
## 1982 - 1980 == 0 6.67336 11.04000 0.604 0.547848
## 1983 - 1980 == 0 -29.18157 12.15792 -2.400 0.019560 *
## 1984 - 1980 == 0 -19.90193 11.39114 -1.747 0.085816 .
## 1985 - 1980 == 0 -9.27748 13.57665 -0.683 0.497066
## 1986 - 1980 == 0 -19.35036 11.39114 -1.699 0.094640 .
## 1987 - 1980 == 0 -9.59250 10.28678 -0.933 0.354875
## 1988 - 1980 == 0 -7.83468 11.45073 -0.684 0.496522
## 1989 - 1980 == 0 -12.15059 10.57357 -1.149 0.255132
## 1990 - 1980 == 0 -13.30936 17.14223 -0.776 0.440610
## 1991 - 1980 == 0 1.73115 12.71073 0.136 0.892130
## 1992 - 1980 == 0 -2.09574 10.51193 -0.199 0.842661
## 1993 - 1980 == 0 -8.72528 11.12874 -0.784 0.436158
## 1994 - 1980 == 0 0.38505 10.88700 0.035 0.971906
## 1995 - 1980 == 0 22.63211 11.42620 1.981 0.052291 .
## 1996 - 1980 == 0 5.70753 11.42620 0.500 0.619276
## 1997 - 1980 == 0 -2.54241 10.64946 -0.239 0.812138
## 1998 - 1980 == 0 23.67982 17.29679 1.369 0.176179
## 1982 - 1981 == 0 -1.60228 16.43022 -0.098 0.922644
## 1983 - 1981 == 0 -37.45721 17.29679 -2.166 0.034399 *
## 1984 - 1981 == 0 -28.17757 16.66532 -1.691 0.096154 .
## 1985 - 1981 == 0 -17.55313 18.12710 -0.968 0.336830
## 1986 - 1981 == 0 -27.62601 16.66532 -1.658 0.102687
## 1987 - 1981 == 0 -17.86814 15.95407 -1.120 0.267263
## 1988 - 1981 == 0 -16.11033 16.79651 -0.959 0.341399
## 1989 - 1981 == 0 -20.42624 16.16342 -1.264 0.211297
## 1990 - 1981 == 0 -21.58500 20.93138 -1.031 0.306643
## 1991 - 1981 == 0 -6.54449 17.69601 -0.370 0.712834
## 1992 - 1981 == 0 -10.37138 16.10151 -0.644 0.521991
## 1993 - 1981 == 0 -17.00092 16.53313 -1.028 0.308009
## 1994 - 1981 == 0 -7.89059 16.37629 -0.482 0.631709
## 1995 - 1981 == 0 14.35647 16.81692 0.854 0.396729
## 1996 - 1981 == 0 -2.56812 16.81692 -0.153 0.879148
## 1997 - 1981 == 0 -10.81805 16.22072 -0.667 0.507418
## 1998 - 1981 == 0 15.40418 21.30941 0.723 0.472609
## 1983 - 1982 == 0 -35.85493 11.10858 -3.228 0.002039 **
## 1984 - 1982 == 0 -26.57529 9.96432 -2.667 0.009861 **
## 1985 - 1982 == 0 -15.95085 12.66577 -1.259 0.212857
## 1986 - 1982 == 0 -26.02373 9.96432 -2.612 0.011408 *
## 1987 - 1982 == 0 -16.26586 8.81152 -1.846 0.069914 .
## 1988 - 1982 == 0 -14.50805 10.30126 -1.408 0.164268
## 1989 - 1982 == 0 -18.82396 9.10569 -2.067 0.043104 *
## 1990 - 1982 == 0 -19.98272 16.43022 -1.216 0.228746
## 1991 - 1982 == 0 -4.94221 11.64016 -0.425 0.672685
## 1992 - 1982 == 0 -8.76910 9.26628 -0.946 0.347832
## 1993 - 1982 == 0 -15.39864 9.67190 -1.592 0.116706
## 1994 - 1982 == 0 -6.28831 9.68212 -0.649 0.518551
## 1995 - 1982 == 0 15.95875 10.29949 1.549 0.126617
## 1996 - 1982 == 0 -0.96584 10.29949 -0.094 0.925606
## 1997 - 1982 == 0 -9.21577 9.39664 -0.981 0.330720
## 1998 - 1982 == 0 17.00646 16.56830 1.026 0.308872
## 1984 - 1983 == 0 9.27964 11.46179 0.810 0.421414
## 1985 - 1983 == 0 19.90409 13.77128 1.445 0.153655
## 1986 - 1983 == 0 9.83121 11.46179 0.858 0.394508
## 1987 - 1983 == 0 19.58907 10.32922 1.896 0.062797 .
## 1988 - 1983 == 0 21.34689 11.38881 1.874 0.065832 .
## 1989 - 1983 == 0 17.03098 10.57989 1.610 0.112791
## 1990 - 1983 == 0 15.87221 17.29679 0.918 0.362544
## 1991 - 1983 == 0 30.91272 12.63195 2.447 0.017394 *
## 1992 - 1983 == 0 27.08583 10.55141 2.567 0.012813 *
## 1993 - 1983 == 0 20.45629 11.13294 1.837 0.071180 .
## 1994 - 1983 == 0 29.56662 10.88394 2.717 0.008643 **
## 1995 - 1983 == 0 51.81368 11.30914 4.582 2.44e-05 ***
## 1996 - 1983 == 0 34.88910 11.30914 3.085 0.003096 **
## 1997 - 1983 == 0 26.63916 10.64424 2.503 0.015116 *
## 1998 - 1983 == 0 52.86139 17.14223 3.084 0.003108 **
## 1985 - 1984 == 0 10.62444 12.96928 0.819 0.415969
## 1986 - 1984 == 0 0.55156 10.46569 0.053 0.958147
## 1987 - 1984 == 0 10.30942 9.32873 1.105 0.273591
## 1988 - 1984 == 0 12.06724 10.69153 1.129 0.263604
## 1989 - 1984 == 0 7.75133 9.62096 0.806 0.423668
## 1990 - 1984 == 0 6.59257 16.66532 0.396 0.693837
## 1991 - 1984 == 0 21.63308 12.00358 1.802 0.076620 .
## 1992 - 1984 == 0 17.80619 9.69120 1.837 0.071195 .
## 1993 - 1984 == 0 11.17665 10.18459 1.097 0.276922
## 1994 - 1984 == 0 20.28698 10.09227 2.010 0.048994 *
## 1995 - 1984 == 0 42.53404 10.68062 3.982 0.000190 ***
## 1996 - 1984 == 0 25.60945 10.68062 2.398 0.019679 *
## 1997 - 1984 == 0 17.35952 9.82504 1.767 0.082423 .
## 1998 - 1984 == 0 43.58175 16.81001 2.593 0.011990 *
## 1986 - 1985 == 0 -10.07288 12.96928 -0.777 0.440455
## 1987 - 1985 == 0 -0.31502 12.04166 -0.026 0.979217
## 1988 - 1985 == 0 1.44280 13.13743 0.110 0.912922
## 1989 - 1985 == 0 -2.87311 12.31769 -0.233 0.816374
## 1990 - 1985 == 0 -4.03187 18.12710 -0.222 0.824753
## 1991 - 1985 == 0 11.00863 14.26948 0.771 0.443501
## 1992 - 1985 == 0 7.18174 12.23634 0.587 0.559497
## 1993 - 1985 == 0 0.55221 12.79898 0.043 0.965732
## 1994 - 1985 == 0 9.66253 12.59573 0.767 0.446065
## 1995 - 1985 == 0 31.90959 13.16352 2.424 0.018430 *
## 1996 - 1985 == 0 14.98501 13.16352 1.138 0.259566
## 1997 - 1985 == 0 6.73508 12.39278 0.543 0.588856
## 1998 - 1985 == 0 32.95730 18.56234 1.775 0.080974 .
## 1987 - 1986 == 0 9.75786 9.32873 1.046 0.299827
## 1988 - 1986 == 0 11.51568 10.69153 1.077 0.285827
## 1989 - 1986 == 0 7.19977 9.62096 0.748 0.457227
## 1990 - 1986 == 0 6.04101 16.66532 0.362 0.718281
## 1991 - 1986 == 0 21.08151 12.00358 1.756 0.084232 .
## 1992 - 1986 == 0 17.25462 9.69120 1.780 0.080153 .
## 1993 - 1986 == 0 10.62509 10.18459 1.043 0.301088
## 1994 - 1986 == 0 19.73541 10.09227 1.955 0.055266 .
## 1995 - 1986 == 0 41.98247 10.68062 3.931 0.000225 ***
## 1996 - 1986 == 0 25.05789 10.68062 2.346 0.022350 *
## 1997 - 1986 == 0 16.80796 9.82504 1.711 0.092385 .
## 1998 - 1986 == 0 43.03019 16.81001 2.560 0.013055 *
## 1988 - 1987 == 0 1.75782 9.36462 0.188 0.851749
## 1989 - 1987 == 0 -2.55809 8.24140 -0.310 0.757354
## 1990 - 1987 == 0 -3.71686 15.95407 -0.233 0.816589
## 1991 - 1987 == 0 11.32365 10.67732 1.061 0.293225
## 1992 - 1987 == 0 7.49676 8.27310 0.906 0.368536
## 1993 - 1987 == 0 0.86722 8.79756 0.099 0.921809
## 1994 - 1987 == 0 9.97755 8.71342 1.145 0.256800
## 1995 - 1987 == 0 32.22461 9.44515 3.412 0.001171 **
## 1996 - 1987 == 0 15.30003 9.44515 1.620 0.110591
## 1997 - 1987 == 0 7.05009 8.33792 0.846 0.401222
## 1998 - 1987 == 0 33.27232 16.03612 2.075 0.042372 *
## 1989 - 1988 == 0 -4.31591 9.62361 -0.448 0.655456
## 1990 - 1988 == 0 -5.47467 16.79651 -0.326 0.745622
## 1991 - 1988 == 0 9.56583 11.53751 0.829 0.410386
## 1992 - 1988 == 0 5.73894 9.60194 0.598 0.552337
## 1993 - 1988 == 0 -0.89059 10.05120 -0.089 0.929695
## 1994 - 1988 == 0 8.21973 9.94532 0.826 0.411855
## 1995 - 1988 == 0 30.46679 10.56650 2.883 0.005483 **
## 1996 - 1988 == 0 13.54221 10.56650 1.282 0.204991
## 1997 - 1988 == 0 5.29228 9.57737 0.553 0.582637
## 1998 - 1988 == 0 31.51450 16.66974 1.891 0.063602 .
## 1990 - 1989 == 0 -1.15876 16.16342 -0.072 0.943091
## 1991 - 1989 == 0 13.88175 10.86908 1.277 0.206543
## 1992 - 1989 == 0 10.05485 8.60023 1.169 0.247050
## 1993 - 1989 == 0 3.42532 9.05878 0.378 0.706697
## 1994 - 1989 == 0 12.53564 9.01126 1.391 0.169416
## 1995 - 1989 == 0 34.78271 9.70911 3.582 0.000690 ***
## 1996 - 1989 == 0 17.85812 9.70911 1.839 0.070902 .
## 1997 - 1989 == 0 9.60819 8.63586 1.113 0.270396
## 1998 - 1989 == 0 35.83042 16.17582 2.215 0.030629 *
## 1991 - 1990 == 0 15.04051 17.69601 0.850 0.398796
## 1992 - 1990 == 0 11.21362 16.10151 0.696 0.488893
## 1993 - 1990 == 0 4.58408 16.53313 0.277 0.782545
## 1994 - 1990 == 0 13.69441 16.37629 0.836 0.406398
## 1995 - 1990 == 0 35.94147 16.81692 2.137 0.036736 *
## 1996 - 1990 == 0 19.01688 16.81692 1.131 0.262707
## 1997 - 1990 == 0 10.76695 16.22072 0.664 0.509418
## 1998 - 1990 == 0 36.98918 21.30941 1.736 0.087817 .
## 1992 - 1991 == 0 -3.82689 10.87482 -0.352 0.726164
## 1993 - 1991 == 0 -10.45643 10.97923 -0.952 0.344788
## 1994 - 1991 == 0 -1.34610 11.13470 -0.121 0.904187
## 1995 - 1991 == 0 20.90096 11.88976 1.758 0.083952 .
## 1996 - 1991 == 0 3.97638 11.88976 0.334 0.739236
## 1997 - 1991 == 0 -4.27356 10.64946 -0.401 0.689654
## 1998 - 1991 == 0 21.94867 17.52595 1.252 0.215383
## 1993 - 1992 == 0 -6.62954 9.14341 -0.725 0.471280
## 1994 - 1992 == 0 2.48079 8.97050 0.277 0.783092
## 1995 - 1992 == 0 24.72785 9.68708 2.553 0.013298 *
## 1996 - 1992 == 0 7.80327 9.68708 0.806 0.423747
## 1997 - 1992 == 0 -0.44667 8.60178 -0.052 0.958762
## 1998 - 1992 == 0 25.77556 16.17879 1.593 0.116466
## 1994 - 1993 == 0 9.11033 9.50529 0.958 0.341748
## 1995 - 1993 == 0 31.35739 10.30850 3.042 0.003505 **
## 1996 - 1993 == 0 14.43280 10.30850 1.400 0.166724
## 1997 - 1993 == 0 6.18287 9.04464 0.684 0.496907
## 1998 - 1993 == 0 32.40510 16.54163 1.959 0.054844 .
## 1995 - 1994 == 0 22.24706 10.03707 2.216 0.030526 *
## 1996 - 1994 == 0 5.32248 10.03707 0.530 0.597907
## 1997 - 1994 == 0 -2.92746 8.99241 -0.326 0.745919
## 1998 - 1994 == 0 23.29477 16.37018 1.423 0.160002
## 1996 - 1995 == 0 -16.92458 10.46569 -1.617 0.111181
## 1997 - 1995 == 0 -25.17452 9.77607 -2.575 0.012548 *
## 1998 - 1995 == 0 1.04771 16.57783 0.063 0.949821
## 1997 - 1996 == 0 -8.24993 9.77607 -0.844 0.402139
## 1998 - 1996 == 0 17.97229 16.57783 1.084 0.282724
## 1998 - 1997 == 0 26.22223 16.21043 1.618 0.111080
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## (Univariate p values reported)
comp_years_ENSO.73_97 <- lmer(Coral_cover ~ Year_F + Thermal_regime + (1|Region), data = ENSO.73_97)
m.comp_ENSO.73_97 <- glht(comp_years_ENSO.73_97, linfct = mcp(Year_F = "Tukey"))
summary(m.comp_ENSO.73_97, test = univariate())
##
## Simultaneous Tests for General Linear Hypotheses
##
## Multiple Comparisons of Means: Tukey Contrasts
##
##
## Fit: lmer(formula = Coral_cover ~ Year_F + Thermal_regime + (1 | Region),
## data = ENSO.73_97)
##
## Linear Hypotheses:
## Estimate Std. Error z value Pr(>|z|)
## 1975 - 1974 == 0 -3.445185 11.191871 -0.308 0.758212
## 1976 - 1974 == 0 10.213685 12.601336 0.811 0.417639
## 1977 - 1974 == 0 0.772940 16.081982 0.048 0.961666
## 1978 - 1974 == 0 -25.039865 13.004321 -1.926 0.054166 .
## 1979 - 1974 == 0 7.674889 12.886959 0.596 0.551473
## 1980 - 1974 == 0 -1.512642 11.560120 -0.131 0.895894
## 1981 - 1974 == 0 9.457940 16.081982 0.588 0.556460
## 1982 - 1974 == 0 10.486135 10.223110 1.026 0.305020
## 1983 - 1974 == 0 -28.692716 11.429261 -2.510 0.012057 *
## 1984 - 1974 == 0 -18.055954 10.535662 -1.714 0.086567 .
## 1985 - 1974 == 0 -8.919951 12.863309 -0.693 0.488033
## 1986 - 1974 == 0 -17.504391 10.535662 -1.661 0.096625 .
## 1987 - 1974 == 0 -8.721286 9.678237 -0.901 0.367523
## 1988 - 1974 == 0 -7.973247 10.762219 -0.741 0.458781
## 1989 - 1974 == 0 -12.132248 9.849498 -1.232 0.218038
## 1990 - 1974 == 0 -12.127060 16.081982 -0.754 0.450803
## 1991 - 1974 == 0 5.151372 11.990928 0.430 0.667482
## 1992 - 1974 == 0 -1.267919 10.045475 -0.126 0.899559
## 1993 - 1974 == 0 -5.459203 10.270773 -0.532 0.595053
## 1994 - 1974 == 0 2.949027 10.389997 0.284 0.776538
## 1995 - 1974 == 0 22.015629 10.809924 2.037 0.041689 *
## 1996 - 1974 == 0 4.678662 10.899444 0.429 0.667736
## 1997 - 1974 == 0 -1.128683 10.254176 -0.110 0.912353
## 1998 - 1974 == 0 18.426838 16.621466 1.109 0.267595
## 1976 - 1975 == 0 13.658870 12.601336 1.084 0.278399
## 1977 - 1975 == 0 4.218126 16.081982 0.262 0.793099
## 1978 - 1975 == 0 -21.594680 13.004321 -1.661 0.096798 .
## 1979 - 1975 == 0 11.120074 12.886959 0.863 0.388196
## 1980 - 1975 == 0 1.932543 11.560120 0.167 0.867234
## 1981 - 1975 == 0 12.903126 16.081982 0.802 0.422360
## 1982 - 1975 == 0 13.931321 10.223110 1.363 0.172968
## 1983 - 1975 == 0 -25.247530 11.429261 -2.209 0.027173 *
## 1984 - 1975 == 0 -14.610769 10.535662 -1.387 0.165505
## 1985 - 1975 == 0 -5.474766 12.863309 -0.426 0.670391
## 1986 - 1975 == 0 -14.059206 10.535662 -1.334 0.182060
## 1987 - 1975 == 0 -5.276100 9.678237 -0.545 0.585650
## 1988 - 1975 == 0 -4.528062 10.762219 -0.421 0.673947
## 1989 - 1975 == 0 -8.687063 9.849498 -0.882 0.377788
## 1990 - 1975 == 0 -8.681874 16.081982 -0.540 0.589300
## 1991 - 1975 == 0 8.596557 11.990928 0.717 0.473422
## 1992 - 1975 == 0 2.177266 10.045475 0.217 0.828410
## 1993 - 1975 == 0 -2.014018 10.270773 -0.196 0.844538
## 1994 - 1975 == 0 6.394212 10.389997 0.615 0.538277
## 1995 - 1975 == 0 25.460814 10.809924 2.355 0.018507 *
## 1996 - 1975 == 0 8.123848 10.899444 0.745 0.456063
## 1997 - 1975 == 0 2.316502 10.254176 0.226 0.821273
## 1998 - 1975 == 0 21.872023 16.621466 1.316 0.188211
## 1977 - 1976 == 0 -9.440744 17.030957 -0.554 0.579354
## 1978 - 1976 == 0 -35.253550 14.262050 -2.472 0.013442 *
## 1979 - 1976 == 0 -2.538796 14.141058 -0.180 0.857519
## 1980 - 1976 == 0 -11.726327 12.948822 -0.906 0.365153
## 1981 - 1976 == 0 -0.755744 17.030957 -0.044 0.964606
## 1982 - 1976 == 0 0.272451 11.718637 0.023 0.981451
## 1983 - 1976 == 0 -38.906400 12.936978 -3.007 0.002635 **
## 1984 - 1976 == 0 -28.269639 12.030409 -2.350 0.018781 *
## 1985 - 1976 == 0 -19.133636 14.083110 -1.359 0.174266
## 1986 - 1976 == 0 -27.718076 12.030409 -2.304 0.021223 *
## 1987 - 1976 == 0 -18.934970 11.319418 -1.673 0.094369 .
## 1988 - 1976 == 0 -18.186932 12.342202 -1.474 0.140601
## 1989 - 1976 == 0 -22.345933 11.468598 -1.948 0.051362 .
## 1990 - 1976 == 0 -22.340744 17.030957 -1.312 0.189597
## 1991 - 1976 == 0 -5.062313 13.495977 -0.375 0.707588
## 1992 - 1976 == 0 -11.481603 11.689973 -0.982 0.326013
## 1993 - 1976 == 0 -15.672888 11.828134 -1.325 0.185154
## 1994 - 1976 == 0 -7.264658 12.006248 -0.605 0.545130
## 1995 - 1976 == 0 11.801944 12.397520 0.952 0.341117
## 1996 - 1976 == 0 -5.535022 12.504596 -0.443 0.658027
## 1997 - 1976 == 0 -11.342368 11.895371 -0.954 0.340331
## 1998 - 1976 == 0 8.213153 17.681669 0.465 0.642289
## 1978 - 1977 == 0 -25.812805 17.101655 -1.509 0.131203
## 1979 - 1977 == 0 6.901949 17.007789 0.406 0.684881
## 1980 - 1977 == 0 -2.285582 16.078972 -0.142 0.886964
## 1981 - 1977 == 0 8.685000 19.384888 0.448 0.654132
## 1982 - 1977 == 0 9.713195 15.439010 0.629 0.529262
## 1983 - 1977 == 0 -29.465656 16.081700 -1.832 0.066915 .
## 1984 - 1977 == 0 -18.828894 15.591871 -1.208 0.227197
## 1985 - 1977 == 0 -9.692892 16.976185 -0.571 0.568020
## 1986 - 1977 == 0 -18.277332 15.591871 -1.172 0.241103
## 1987 - 1977 == 0 -9.494226 14.975579 -0.634 0.526094
## 1988 - 1977 == 0 -8.746188 15.632770 -0.559 0.575836
## 1989 - 1977 == 0 -12.905188 15.110662 -0.854 0.393080
## 1990 - 1977 == 0 -12.900000 19.384888 -0.665 0.505752
## 1991 - 1977 == 0 4.378432 16.635386 0.263 0.792397
## 1992 - 1977 == 0 -2.040859 15.134351 -0.135 0.892731
## 1993 - 1977 == 0 -6.232144 15.414473 -0.404 0.685989
## 1994 - 1977 == 0 2.176086 15.379144 0.141 0.887478
## 1995 - 1977 == 0 21.242689 15.660154 1.356 0.174946
## 1996 - 1977 == 0 3.905722 15.810600 0.247 0.804884
## 1997 - 1977 == 0 -1.901623 15.284004 -0.124 0.900983
## 1998 - 1977 == 0 17.653898 20.137322 0.877 0.380663
## 1979 - 1978 == 0 32.714754 14.116619 2.317 0.020478 *
## 1980 - 1978 == 0 23.527223 12.963119 1.815 0.069534 .
## 1981 - 1978 == 0 34.497805 17.101655 2.017 0.043673 *
## 1982 - 1978 == 0 35.526000 12.180022 2.917 0.003537 **
## 1983 - 1978 == 0 -3.652851 12.918112 -0.283 0.777353
## 1984 - 1978 == 0 6.983911 12.373073 0.564 0.572452
## 1985 - 1978 == 0 16.119914 14.169928 1.138 0.255282
## 1986 - 1978 == 0 7.535474 12.373073 0.609 0.542510
## 1987 - 1978 == 0 16.318579 11.527100 1.416 0.156872
## 1988 - 1978 == 0 17.066618 12.341141 1.383 0.166694
## 1989 - 1978 == 0 12.907617 11.691774 1.104 0.269597
## 1990 - 1978 == 0 12.912805 17.101655 0.755 0.450212
## 1991 - 1978 == 0 30.191237 13.491233 2.238 0.025231 *
## 1992 - 1978 == 0 23.771946 11.722637 2.028 0.042574 *
## 1993 - 1978 == 0 19.580662 12.102902 1.618 0.105695
## 1994 - 1978 == 0 27.988892 11.770175 2.378 0.017409 *
## 1995 - 1978 == 0 47.055494 12.040342 3.908 9.30e-05 ***
## 1996 - 1978 == 0 29.718527 12.151383 2.446 0.014457 *
## 1997 - 1978 == 0 23.911182 11.670049 2.049 0.040468 *
## 1998 - 1978 == 0 43.466703 17.574047 2.473 0.013385 *
## 1980 - 1979 == 0 -9.187531 12.620705 -0.728 0.466630
## 1981 - 1979 == 0 1.783051 17.007789 0.105 0.916505
## 1982 - 1979 == 0 2.811246 12.047057 0.233 0.815485
## 1983 - 1979 == 0 -36.367605 12.591081 -2.888 0.003873 **
## 1984 - 1979 == 0 -25.730843 12.244988 -2.101 0.035611 *
## 1985 - 1979 == 0 -16.594841 14.051249 -1.181 0.237594
## 1986 - 1979 == 0 -25.179281 12.244988 -2.056 0.039754 *
## 1987 - 1979 == 0 -16.396175 11.278464 -1.454 0.146013
## 1988 - 1979 == 0 -15.648137 12.051288 -1.298 0.194129
## 1989 - 1979 == 0 -19.807137 11.436319 -1.732 0.083282 .
## 1990 - 1979 == 0 -19.801949 17.007789 -1.164 0.244308
## 1991 - 1979 == 0 -2.523517 13.419856 -0.188 0.850843
## 1992 - 1979 == 0 -8.942808 11.463717 -0.780 0.435334
## 1993 - 1979 == 0 -13.134092 11.991057 -1.095 0.273375
## 1994 - 1979 == 0 -4.725862 11.904834 -0.397 0.691390
## 1995 - 1979 == 0 14.340740 12.050292 1.190 0.234017
## 1996 - 1979 == 0 -2.996227 12.158209 -0.246 0.805344
## 1997 - 1979 == 0 -8.803572 11.780457 -0.747 0.454881
## 1998 - 1979 == 0 10.751949 17.007789 0.632 0.527271
## 1981 - 1980 == 0 10.970582 16.078972 0.682 0.495053
## 1982 - 1980 == 0 11.998777 10.478966 1.145 0.252195
## 1983 - 1980 == 0 -27.180074 11.338424 -2.397 0.016522 *
## 1984 - 1980 == 0 -16.543312 10.695603 -1.547 0.121926
## 1985 - 1980 == 0 -7.407309 12.634881 -0.586 0.557702
## 1986 - 1980 == 0 -15.991749 10.695603 -1.495 0.134870
## 1987 - 1980 == 0 -7.208644 9.659098 -0.746 0.455483
## 1988 - 1980 == 0 -6.460605 10.699788 -0.604 0.545972
## 1989 - 1980 == 0 -10.619606 9.857462 -1.077 0.281339
## 1990 - 1980 == 0 -10.614418 16.078972 -0.660 0.509162
## 1991 - 1980 == 0 6.664014 12.137682 0.549 0.582981
## 1992 - 1980 == 0 0.244723 9.874423 0.025 0.980228
## 1993 - 1980 == 0 -3.946561 10.535916 -0.375 0.707972
## 1994 - 1980 == 0 4.461669 10.344569 0.431 0.666246
## 1995 - 1980 == 0 23.528271 10.714130 2.196 0.028092 *
## 1996 - 1980 == 0 6.191304 10.671091 0.580 0.561784
## 1997 - 1980 == 0 0.383959 10.295711 0.037 0.970251
## 1998 - 1980 == 0 19.939480 16.210705 1.230 0.218690
## 1982 - 1981 == 0 1.028195 15.439010 0.067 0.946902
## 1983 - 1981 == 0 -38.150656 16.081700 -2.372 0.017678 *
## 1984 - 1981 == 0 -27.513894 15.591871 -1.765 0.077626 .
## 1985 - 1981 == 0 -18.377892 16.976185 -1.083 0.279000
## 1986 - 1981 == 0 -26.962332 15.591871 -1.729 0.083763 .
## 1987 - 1981 == 0 -18.179226 14.975579 -1.214 0.224776
## 1988 - 1981 == 0 -17.431188 15.632770 -1.115 0.264833
## 1989 - 1981 == 0 -21.590188 15.110662 -1.429 0.153060
## 1990 - 1981 == 0 -21.585000 19.384888 -1.113 0.265495
## 1991 - 1981 == 0 -4.306568 16.635386 -0.259 0.795728
## 1992 - 1981 == 0 -10.725859 15.134351 -0.709 0.478505
## 1993 - 1981 == 0 -14.917144 15.414473 -0.968 0.333176
## 1994 - 1981 == 0 -6.508914 15.379144 -0.423 0.672128
## 1995 - 1981 == 0 12.557689 15.660154 0.802 0.422618
## 1996 - 1981 == 0 -4.779278 15.810600 -0.302 0.762436
## 1997 - 1981 == 0 -10.586623 15.284004 -0.693 0.488523
## 1998 - 1981 == 0 8.968898 20.137322 0.445 0.656040
## 1983 - 1982 == 0 -39.178851 10.499147 -3.732 0.000190 ***
## 1984 - 1982 == 0 -28.542089 9.383017 -3.042 0.002351 **
## 1985 - 1982 == 0 -19.406087 11.854611 -1.637 0.101629
## 1986 - 1982 == 0 -27.990527 9.383017 -2.983 0.002853 **
## 1987 - 1982 == 0 -19.207421 8.425735 -2.280 0.022631 *
## 1988 - 1982 == 0 -18.459383 9.752771 -1.893 0.058394 .
## 1989 - 1982 == 0 -22.618383 8.631091 -2.621 0.008778 **
## 1990 - 1982 == 0 -22.613195 15.439010 -1.465 0.143008
## 1991 - 1982 == 0 -5.334763 11.097988 -0.481 0.630732
## 1992 - 1982 == 0 -11.754054 8.845054 -1.329 0.183886
## 1993 - 1982 == 0 -15.945339 9.185269 -1.736 0.082569 .
## 1994 - 1982 == 0 -7.537109 9.247990 -0.815 0.415072
## 1995 - 1982 == 0 11.529494 9.810432 1.175 0.239904
## 1996 - 1982 == 0 -5.807473 9.790562 -0.593 0.553067
## 1997 - 1982 == 0 -11.614818 9.165586 -1.267 0.205076
## 1998 - 1982 == 0 7.940702 15.957946 0.498 0.618765
## 1984 - 1983 == 0 10.636762 10.723316 0.992 0.321232
## 1985 - 1983 == 0 19.772764 12.860163 1.538 0.124166
## 1986 - 1983 == 0 11.188324 10.723316 1.043 0.296780
## 1987 - 1983 == 0 19.971430 9.673054 2.065 0.038957 *
## 1988 - 1983 == 0 20.719468 10.564616 1.961 0.049854 *
## 1989 - 1983 == 0 16.560468 9.850705 1.681 0.092735 .
## 1990 - 1983 == 0 16.565656 16.081700 1.030 0.302966
## 1991 - 1983 == 0 33.844088 11.955465 2.831 0.004642 **
## 1992 - 1983 == 0 27.424797 9.896252 2.771 0.005584 **
## 1993 - 1983 == 0 23.233512 10.420777 2.230 0.025778 *
## 1994 - 1983 == 0 31.641742 10.336929 3.061 0.002206 **
## 1995 - 1983 == 0 50.708345 10.569381 4.798 1.61e-06 ***
## 1996 - 1983 == 0 33.371378 10.659236 3.131 0.001744 **
## 1997 - 1983 == 0 27.564033 10.190731 2.705 0.006834 **
## 1998 - 1983 == 0 47.119554 16.161852 2.915 0.003551 **
## 1985 - 1984 == 0 9.136002 12.040605 0.759 0.447993
## 1986 - 1984 == 0 0.551563 9.692444 0.057 0.954620
## 1987 - 1984 == 0 9.334668 8.745813 1.067 0.285823
## 1988 - 1984 == 0 10.082706 9.999131 1.008 0.313283
## 1989 - 1984 == 0 5.923706 8.950319 0.662 0.508072
## 1990 - 1984 == 0 5.928894 15.591871 0.380 0.703756
## 1991 - 1984 == 0 23.207326 11.320304 2.050 0.040358 *
## 1992 - 1984 == 0 16.788035 9.108772 1.843 0.065320 .
## 1993 - 1984 == 0 12.596751 9.504960 1.325 0.185078
## 1994 - 1984 == 0 21.004981 9.498702 2.211 0.027011 *
## 1995 - 1984 == 0 40.071583 10.049715 3.987 6.68e-05 ***
## 1996 - 1984 == 0 22.734616 10.022864 2.268 0.023312 *
## 1997 - 1984 == 0 16.927271 9.425446 1.796 0.072509 .
## 1998 - 1984 == 0 36.482792 16.109900 2.265 0.023536 *
## 1986 - 1985 == 0 -8.584440 12.040605 -0.713 0.475872
## 1987 - 1985 == 0 0.198666 11.253724 0.018 0.985915
## 1988 - 1985 == 0 0.946704 12.266018 0.077 0.938480
## 1989 - 1985 == 0 -3.212296 11.441145 -0.281 0.778889
## 1990 - 1985 == 0 -3.207108 16.976185 -0.189 0.850157
## 1991 - 1985 == 0 14.071324 13.454712 1.046 0.295640
## 1992 - 1985 == 0 7.652033 11.451315 0.668 0.503991
## 1993 - 1985 == 0 3.460748 11.970573 0.289 0.772501
## 1994 - 1985 == 0 11.868978 11.782101 1.007 0.313755
## 1995 - 1985 == 0 30.935580 12.303185 2.514 0.011922 *
## 1996 - 1985 == 0 13.598614 12.236934 1.111 0.266449
## 1997 - 1985 == 0 7.791269 11.781024 0.661 0.508394
## 1998 - 1985 == 0 27.346789 17.590912 1.555 0.120042
## 1987 - 1986 == 0 8.783106 8.745813 1.004 0.315251
## 1988 - 1986 == 0 9.531144 9.999131 0.953 0.340490
## 1989 - 1986 == 0 5.372143 8.950319 0.600 0.548361
## 1990 - 1986 == 0 5.377332 15.591871 0.345 0.730184
## 1991 - 1986 == 0 22.655763 11.320304 2.001 0.045356 *
## 1992 - 1986 == 0 16.236473 9.108772 1.783 0.074666 .
## 1993 - 1986 == 0 12.045188 9.504960 1.267 0.205065
## 1994 - 1986 == 0 20.453418 9.498702 2.153 0.031296 *
## 1995 - 1986 == 0 39.520020 10.049715 3.932 8.41e-05 ***
## 1996 - 1986 == 0 22.183054 10.022864 2.213 0.026881 *
## 1997 - 1986 == 0 16.375709 9.425446 1.737 0.082318 .
## 1998 - 1986 == 0 35.931229 16.109900 2.230 0.025722 *
## 1988 - 1987 == 0 0.748038 8.755306 0.085 0.931913
## 1989 - 1987 == 0 -3.410962 7.693660 -0.443 0.657515
## 1990 - 1987 == 0 -3.405774 14.975579 -0.227 0.820096
## 1991 - 1987 == 0 13.872658 10.192240 1.361 0.173482
## 1992 - 1987 == 0 7.453367 7.665205 0.972 0.330870
## 1993 - 1987 == 0 3.262082 8.328932 0.392 0.695312
## 1994 - 1987 == 0 11.670313 8.369579 1.394 0.163205
## 1995 - 1987 == 0 30.736915 8.918004 3.447 0.000568 ***
## 1996 - 1987 == 0 13.399948 8.902028 1.505 0.132255
## 1997 - 1987 == 0 7.592603 8.135978 0.933 0.350710
## 1998 - 1987 == 0 27.148124 15.260910 1.779 0.075251 .
## 1989 - 1988 == 0 -4.159000 8.945157 -0.465 0.641971
## 1990 - 1988 == 0 -4.153812 15.632770 -0.266 0.790461
## 1991 - 1988 == 0 13.124620 10.966747 1.197 0.231398
## 1992 - 1988 == 0 6.705329 8.991292 0.746 0.455814
## 1993 - 1988 == 0 2.514044 9.453105 0.266 0.790278
## 1994 - 1988 == 0 10.922274 9.519842 1.147 0.251251
## 1995 - 1988 == 0 29.988876 9.886816 3.033 0.002420 **
## 1996 - 1988 == 0 12.651910 9.965497 1.270 0.204237
## 1997 - 1988 == 0 6.844565 9.229244 0.742 0.458319
## 1998 - 1988 == 0 26.400085 15.774645 1.674 0.094214 .
## 1990 - 1989 == 0 0.005188 15.110662 0.000 0.999726
## 1991 - 1989 == 0 17.283620 10.339771 1.672 0.094610 .
## 1992 - 1989 == 0 10.864329 8.045393 1.350 0.176894
## 1993 - 1989 == 0 6.673045 8.525934 0.783 0.433817
## 1994 - 1989 == 0 15.081275 8.591925 1.755 0.079211 .
## 1995 - 1989 == 0 34.147877 9.111962 3.748 0.000179 ***
## 1996 - 1989 == 0 16.810910 9.101530 1.847 0.064741 .
## 1997 - 1989 == 0 11.003565 8.382689 1.313 0.189300
## 1998 - 1989 == 0 30.559086 15.362374 1.989 0.046677 *
## 1991 - 1990 == 0 17.278432 16.635386 1.039 0.298965
## 1992 - 1990 == 0 10.859141 15.134351 0.718 0.473056
## 1993 - 1990 == 0 6.667856 15.414473 0.433 0.665326
## 1994 - 1990 == 0 15.076086 15.379144 0.980 0.326941
## 1995 - 1990 == 0 34.142689 15.660154 2.180 0.029241 *
## 1996 - 1990 == 0 16.805722 15.810600 1.063 0.287809
## 1997 - 1990 == 0 10.998377 15.284004 0.720 0.471771
## 1998 - 1990 == 0 30.553898 20.137322 1.517 0.129197
## 1992 - 1991 == 0 -6.419291 10.395305 -0.618 0.536893
## 1993 - 1991 == 0 -10.610575 10.492563 -1.011 0.311898
## 1994 - 1991 == 0 -2.202345 10.714982 -0.206 0.837151
## 1995 - 1991 == 0 16.864257 11.314032 1.491 0.136077
## 1996 - 1991 == 0 -0.472710 11.343724 -0.042 0.966761
## 1997 - 1991 == 0 -6.280055 9.966742 -0.630 0.528628
## 1998 - 1991 == 0 13.275466 16.917590 0.785 0.432621
## 1993 - 1992 == 0 -4.191285 8.680851 -0.483 0.629224
## 1994 - 1992 == 0 4.216946 8.634387 0.488 0.625274
## 1995 - 1992 == 0 23.283548 9.158855 2.542 0.011016 *
## 1996 - 1992 == 0 5.946581 9.135250 0.651 0.515079
## 1997 - 1992 == 0 0.139236 8.393537 0.017 0.986765
## 1998 - 1992 == 0 19.694757 15.379886 1.281 0.200351
## 1994 - 1993 == 0 8.408230 8.911097 0.944 0.345390
## 1995 - 1993 == 0 27.474832 9.709164 2.830 0.004658 **
## 1996 - 1993 == 0 10.137866 9.779485 1.037 0.299901
## 1997 - 1993 == 0 4.330520 8.797278 0.492 0.622538
## 1998 - 1993 == 0 23.886041 15.897194 1.503 0.132960
## 1995 - 1994 == 0 19.066602 9.447903 2.018 0.043583 *
## 1996 - 1994 == 0 1.729635 9.430744 0.183 0.854481
## 1997 - 1994 == 0 -4.077710 8.740245 -0.467 0.640826
## 1998 - 1994 == 0 15.477811 15.801517 0.980 0.327326
## 1996 - 1995 == 0 -17.336967 9.774128 -1.774 0.076103 .
## 1997 - 1995 == 0 -23.144312 9.304186 -2.488 0.012864 *
## 1998 - 1995 == 0 -3.588791 15.745937 -0.228 0.819710
## 1997 - 1996 == 0 -5.807345 9.364289 -0.620 0.535153
## 1998 - 1996 == 0 13.748176 15.761500 0.872 0.383065
## 1998 - 1997 == 0 19.555521 15.707114 1.245 0.213128
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## (Univariate p values reported)
ENSO.83_97$Year_F <- as.factor(ENSO.83_97$Year)
comp_years_ENSO.83_97 <- aov(Coral_cover ~ Year_F + Thermal_regime, data = ENSO.83_97)
m.comp_ENSO.83_97 <- glht(comp_years_ENSO.83_97, linfct = mcp(Year_F = "Tukey"))
summary(m.comp_ENSO.83_97, test = univariate())
##
## Simultaneous Tests for General Linear Hypotheses
##
## Multiple Comparisons of Means: Tukey Contrasts
##
##
## Fit: aov(formula = Coral_cover ~ Year_F + Thermal_regime, data = ENSO.83_97)
##
## Linear Hypotheses:
## Estimate Std. Error t value Pr(>|t|)
## 1984 - 1983 == 0 12.351578 10.375355 1.190 0.239968
## 1985 - 1983 == 0 21.853429 12.451269 1.755 0.085898 .
## 1986 - 1983 == 0 12.903140 10.375355 1.244 0.219935
## 1987 - 1983 == 0 21.659577 9.319358 2.324 0.024588 *
## 1988 - 1983 == 0 21.654675 10.234189 2.116 0.039796 *
## 1989 - 1983 == 0 19.121678 9.546497 2.003 0.051090 .
## 1990 - 1983 == 0 17.821554 15.603472 1.142 0.259301
## 1991 - 1983 == 0 31.408598 11.354258 2.766 0.008138 **
## 1992 - 1983 == 0 27.940801 9.494612 2.943 0.005082 **
## 1993 - 1983 == 0 23.016432 10.063858 2.287 0.026843 *
## 1994 - 1983 == 0 30.202718 9.785636 3.086 0.003426 **
## 1995 - 1983 == 0 51.570014 10.163303 5.074 6.85e-06 ***
## 1996 - 1983 == 0 34.645431 10.163303 3.409 0.001366 **
## 1997 - 1983 == 0 27.374437 9.570674 2.860 0.006346 **
## 1998 - 1983 == 0 51.886723 15.418676 3.365 0.001551 **
## 1985 - 1984 == 0 9.501851 11.732643 0.810 0.422188
## 1986 - 1984 == 0 0.551563 9.403948 0.059 0.953483
## 1987 - 1984 == 0 9.307999 8.393339 1.109 0.273206
## 1988 - 1984 == 0 9.303098 9.673243 0.962 0.341211
## 1989 - 1984 == 0 6.770100 8.654255 0.782 0.438055
## 1990 - 1984 == 0 5.469976 15.036262 0.364 0.717685
## 1991 - 1984 == 0 19.057020 10.836676 1.759 0.085303 .
## 1992 - 1984 == 0 15.589224 8.769831 1.778 0.082081 .
## 1993 - 1984 == 0 10.664854 9.156204 1.165 0.250118
## 1994 - 1984 == 0 17.851140 9.129581 1.955 0.056638 .
## 1995 - 1984 == 0 39.218436 9.691117 4.047 0.000197 ***
## 1996 - 1984 == 0 22.293853 9.691117 2.300 0.026008 *
## 1997 - 1984 == 0 15.022859 8.885050 1.691 0.097641 .
## 1998 - 1984 == 0 39.535145 15.199628 2.601 0.012458 *
## 1986 - 1985 == 0 -8.950288 11.732643 -0.763 0.449445
## 1987 - 1985 == 0 -0.193851 10.880330 -0.018 0.985862
## 1988 - 1985 == 0 -0.198753 11.878661 -0.017 0.986723
## 1989 - 1985 == 0 -2.731750 11.148223 -0.245 0.807516
## 1990 - 1985 == 0 -4.031875 16.288116 -0.248 0.805596
## 1991 - 1985 == 0 9.555169 12.906122 0.740 0.462844
## 1992 - 1985 == 0 6.087373 11.030341 0.552 0.583705
## 1993 - 1985 == 0 1.163004 11.599831 0.100 0.920573
## 1994 - 1985 == 0 8.349289 11.367353 0.734 0.466374
## 1995 - 1985 == 0 29.716586 11.930022 2.491 0.016413 *
## 1996 - 1985 == 0 12.792002 11.930022 1.072 0.289200
## 1997 - 1985 == 0 5.521009 11.187549 0.493 0.624010
## 1998 - 1985 == 0 30.033295 16.807767 1.787 0.080549 .
## 1987 - 1986 == 0 8.756437 8.393339 1.043 0.302280
## 1988 - 1986 == 0 8.751535 9.673243 0.905 0.370330
## 1989 - 1986 == 0 6.218538 8.654255 0.719 0.476051
## 1990 - 1986 == 0 4.918413 15.036262 0.327 0.745074
## 1991 - 1986 == 0 18.505457 10.836676 1.708 0.094439 .
## 1992 - 1986 == 0 15.037661 8.769831 1.715 0.093130 .
## 1993 - 1986 == 0 10.113292 9.156204 1.105 0.275108
## 1994 - 1986 == 0 17.299577 9.129581 1.895 0.064402 .
## 1995 - 1986 == 0 38.666874 9.691117 3.990 0.000235 ***
## 1996 - 1986 == 0 21.742290 9.691117 2.244 0.029720 *
## 1997 - 1986 == 0 14.471297 8.885050 1.629 0.110202
## 1998 - 1986 == 0 38.983583 15.199628 2.565 0.013653 *
## 1988 - 1987 == 0 -0.004902 8.444989 -0.001 0.999539
## 1989 - 1987 == 0 -2.537899 7.407768 -0.343 0.733460
## 1990 - 1987 == 0 -3.838024 14.381093 -0.267 0.790755
## 1991 - 1987 == 0 9.749021 9.617622 1.014 0.316046
## 1992 - 1987 == 0 6.281224 7.459867 0.842 0.404143
## 1993 - 1987 == 0 1.356855 7.915355 0.171 0.864645
## 1994 - 1987 == 0 8.543141 7.854396 1.088 0.282399
## 1995 - 1987 == 0 29.910437 8.539886 3.502 0.001038 **
## 1996 - 1987 == 0 12.985854 8.539886 1.521 0.135202
## 1997 - 1987 == 0 5.714860 7.513916 0.761 0.450795
## 1998 - 1987 == 0 30.227146 14.473044 2.089 0.042314 *
## 1989 - 1988 == 0 -2.532997 8.679641 -0.292 0.771727
## 1990 - 1988 == 0 -3.833122 15.150473 -0.253 0.801393
## 1991 - 1988 == 0 9.753922 10.369343 0.941 0.351798
## 1992 - 1988 == 0 6.286126 8.639110 0.728 0.470524
## 1993 - 1988 == 0 1.361757 9.084485 0.150 0.881499
## 1994 - 1988 == 0 8.548042 8.940295 0.956 0.344007
## 1995 - 1988 == 0 29.915339 9.498092 3.150 0.002871 **
## 1996 - 1988 == 0 12.990755 9.498092 1.368 0.178044
## 1997 - 1988 == 0 5.719762 8.609504 0.664 0.509780
## 1998 - 1988 == 0 30.232048 14.998026 2.016 0.049691 *
## 1990 - 1989 == 0 -1.300125 14.584825 -0.089 0.929356
## 1991 - 1989 == 0 12.286919 9.787894 1.255 0.215703
## 1992 - 1989 == 0 8.819123 7.765025 1.136 0.261944
## 1993 - 1989 == 0 3.894754 8.143592 0.478 0.634730
## 1994 - 1989 == 0 11.081039 8.129358 1.363 0.179490
## 1995 - 1989 == 0 32.448336 8.775695 3.698 0.000579 ***
## 1996 - 1989 == 0 15.523753 8.775695 1.769 0.083532 .
## 1997 - 1989 == 0 8.252759 7.788304 1.060 0.294843
## 1998 - 1989 == 0 32.765045 14.592673 2.245 0.029597 *
## 1991 - 1990 == 0 13.587044 15.968790 0.851 0.399260
## 1992 - 1990 == 0 10.119248 14.494919 0.698 0.488614
## 1993 - 1990 == 0 5.194879 14.932861 0.348 0.729515
## 1994 - 1990 == 0 12.381164 14.752999 0.839 0.405680
## 1995 - 1990 == 0 33.748461 15.190776 2.222 0.031268 *
## 1996 - 1990 == 0 16.823877 15.190776 1.108 0.273833
## 1997 - 1990 == 0 9.552884 14.614907 0.654 0.516598
## 1998 - 1990 == 0 34.065170 19.259680 1.769 0.083569 .
## 1992 - 1991 == 0 -3.467796 9.790236 -0.354 0.724800
## 1993 - 1991 == 0 -8.392165 9.902026 -0.848 0.401094
## 1994 - 1991 == 0 -1.205880 10.014648 -0.120 0.904681
## 1995 - 1991 == 0 20.161416 10.688016 1.886 0.065569 .
## 1996 - 1991 == 0 3.236833 10.688016 0.303 0.763372
## 1997 - 1991 == 0 -4.034160 9.577556 -0.421 0.675564
## 1998 - 1991 == 0 20.478126 15.763214 1.299 0.200381
## 1993 - 1992 == 0 -4.924369 8.277348 -0.595 0.554812
## 1994 - 1992 == 0 2.261916 8.062383 0.281 0.780312
## 1995 - 1992 == 0 23.629213 8.730875 2.706 0.009511 **
## 1996 - 1992 == 0 6.704629 8.730875 0.768 0.446458
## 1997 - 1992 == 0 -0.566364 7.731752 -0.073 0.941923
## 1998 - 1992 == 0 23.945922 14.586686 1.642 0.107487
## 1994 - 1993 == 0 7.186285 8.596189 0.836 0.407485
## 1995 - 1993 == 0 28.553582 9.335014 3.059 0.003699 **
## 1996 - 1993 == 0 11.628999 9.335014 1.246 0.219168
## 1997 - 1993 == 0 4.358005 8.178387 0.533 0.596690
## 1998 - 1993 == 0 28.870291 14.934120 1.933 0.059383 .
## 1995 - 1994 == 0 21.367296 9.033108 2.365 0.022280 *
## 1996 - 1994 == 0 4.442713 9.033108 0.492 0.625182
## 1997 - 1994 == 0 -2.828280 8.080350 -0.350 0.727922
## 1998 - 1994 == 0 21.684006 14.744958 1.471 0.148206
## 1996 - 1995 == 0 -16.924583 9.403948 -1.800 0.078464 .
## 1997 - 1995 == 0 -24.195577 8.799034 -2.750 0.008496 **
## 1998 - 1995 == 0 0.316709 14.905050 0.021 0.983139
## 1997 - 1996 == 0 -7.270994 8.799034 -0.826 0.412877
## 1998 - 1996 == 0 17.241292 14.905050 1.157 0.253350
## 1998 - 1997 == 0 24.512286 14.601377 1.679 0.099979 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## (Univariate p values reported)
comp_years_ENSO.83_97 <- lmer(Coral_cover ~ Year_F + Thermal_regime + (1|Region), data = ENSO.83_97)
m.comp_ENSO.83_97 <- glht(comp_years_ENSO.83_97, linfct = mcp(Year_F = "Tukey"))
summary(m.comp_ENSO.83_97, test = univariate())
##
## Simultaneous Tests for General Linear Hypotheses
##
## Multiple Comparisons of Means: Tukey Contrasts
##
##
## Fit: lmer(formula = Coral_cover ~ Year_F + Thermal_regime + (1 | Region),
## data = ENSO.83_97)
##
## Linear Hypotheses:
## Estimate Std. Error z value Pr(>|z|)
## 1984 - 1983 == 0 13.5715 9.3127 1.457 0.145033
## 1985 - 1983 == 0 21.3361 11.1455 1.914 0.055579 .
## 1986 - 1983 == 0 14.1230 9.3127 1.517 0.129386
## 1987 - 1983 == 0 21.9876 8.3733 2.626 0.008642 **
## 1988 - 1983 == 0 20.7808 9.0881 2.287 0.022219 *
## 1989 - 1983 == 0 18.4963 8.5201 2.171 0.029939 *
## 1990 - 1983 == 0 18.0313 13.9183 1.296 0.195146
## 1991 - 1983 == 0 33.8098 10.3383 3.270 0.001074 **
## 1992 - 1983 == 0 28.1589 8.5416 3.297 0.000978 ***
## 1993 - 1983 == 0 25.8114 9.0354 2.857 0.004281 **
## 1994 - 1983 == 0 31.5048 8.9542 3.518 0.000434 ***
## 1995 - 1983 == 0 49.3879 9.1119 5.420 5.95e-08 ***
## 1996 - 1983 == 0 32.0998 9.2215 3.481 0.000500 ***
## 1997 - 1983 == 0 26.8883 8.8386 3.042 0.002349 **
## 1998 - 1983 == 0 46.7503 13.9792 3.344 0.000825 ***
## 1985 - 1984 == 0 7.7647 10.4219 0.745 0.456251
## 1986 - 1984 == 0 0.5516 8.3345 0.066 0.947236
## 1987 - 1984 == 0 8.4161 7.5505 1.115 0.265000
## 1988 - 1984 == 0 7.2093 8.6782 0.831 0.406122
## 1989 - 1984 == 0 4.9248 7.7129 0.639 0.523137
## 1990 - 1984 == 0 4.4598 13.5298 0.330 0.741682
## 1991 - 1984 == 0 20.2384 9.8276 2.059 0.039462 *
## 1992 - 1984 == 0 14.5875 7.9054 1.845 0.065001 .
## 1993 - 1984 == 0 12.2399 8.1927 1.494 0.135176
## 1994 - 1984 == 0 17.9333 8.2490 2.174 0.029706 *
## 1995 - 1984 == 0 35.8164 8.7538 4.092 4.28e-05 ***
## 1996 - 1984 == 0 18.5284 8.7273 2.123 0.033752 *
## 1997 - 1984 == 0 13.3169 8.2289 1.618 0.105596
## 1998 - 1984 == 0 33.1788 14.0467 2.362 0.018174 *
## 1986 - 1985 == 0 -7.2131 10.4219 -0.692 0.488868
## 1987 - 1985 == 0 0.6515 9.7578 0.067 0.946769
## 1988 - 1985 == 0 -0.5554 10.6302 -0.052 0.958334
## 1989 - 1985 == 0 -2.8399 9.9144 -0.286 0.774542
## 1990 - 1985 == 0 -3.3049 14.6533 -0.226 0.821561
## 1991 - 1985 == 0 12.4737 11.6918 1.067 0.286028
## 1992 - 1985 == 0 6.8228 9.9104 0.688 0.491170
## 1993 - 1985 == 0 4.4753 10.3892 0.431 0.666643
## 1994 - 1985 == 0 10.1686 10.1949 0.997 0.318556
## 1995 - 1985 == 0 28.0518 10.6865 2.625 0.008665 **
## 1996 - 1985 == 0 10.7637 10.6331 1.012 0.311402
## 1997 - 1985 == 0 5.5522 10.2443 0.542 0.587833
## 1998 - 1985 == 0 25.4142 15.3383 1.657 0.097538 .
## 1987 - 1986 == 0 7.8646 7.5505 1.042 0.297597
## 1988 - 1986 == 0 6.6577 8.6782 0.767 0.442975
## 1989 - 1986 == 0 4.3732 7.7129 0.567 0.570711
## 1990 - 1986 == 0 3.9082 13.5298 0.289 0.772688
## 1991 - 1986 == 0 19.6868 9.8276 2.003 0.045155 *
## 1992 - 1986 == 0 14.0359 7.9054 1.775 0.075819 .
## 1993 - 1986 == 0 11.6884 8.1927 1.427 0.153674
## 1994 - 1986 == 0 17.3817 8.2490 2.107 0.035107 *
## 1995 - 1986 == 0 35.2649 8.7538 4.029 5.61e-05 ***
## 1996 - 1986 == 0 17.9768 8.7273 2.060 0.039415 *
## 1997 - 1986 == 0 12.7653 8.2289 1.551 0.120835
## 1998 - 1986 == 0 32.6273 14.0467 2.323 0.020191 *
## 1988 - 1987 == 0 -1.2068 7.5723 -0.159 0.873372
## 1989 - 1987 == 0 -3.4913 6.6277 -0.527 0.598345
## 1990 - 1987 == 0 -3.9563 13.0035 -0.304 0.760935
## 1991 - 1987 == 0 11.8222 8.8484 1.336 0.181519
## 1992 - 1987 == 0 6.1713 6.6133 0.933 0.350732
## 1993 - 1987 == 0 3.8238 7.2018 0.531 0.595455
## 1994 - 1987 == 0 9.5172 7.2815 1.307 0.191201
## 1995 - 1987 == 0 27.4003 7.7452 3.538 0.000404 ***
## 1996 - 1987 == 0 10.1122 7.7336 1.308 0.191021
## 1997 - 1987 == 0 4.9007 7.1107 0.689 0.490697
## 1998 - 1987 == 0 24.7627 13.2586 1.868 0.061808 .
## 1989 - 1988 == 0 -2.2845 7.7303 -0.296 0.767594
## 1990 - 1988 == 0 -2.7495 13.5351 -0.203 0.839027
## 1991 - 1988 == 0 13.0291 9.4892 1.373 0.169741
## 1992 - 1988 == 0 7.3782 7.7557 0.951 0.341441
## 1993 - 1988 == 0 5.0306 8.2014 0.613 0.539618
## 1994 - 1988 == 0 10.7240 8.2567 1.299 0.194004
## 1995 - 1988 == 0 28.6071 8.5264 3.355 0.000793 ***
## 1996 - 1988 == 0 11.3191 8.6232 1.313 0.189308
## 1997 - 1988 == 0 6.1076 8.0117 0.762 0.445863
## 1998 - 1988 == 0 25.9696 13.6574 1.902 0.057236 .
## 1990 - 1989 == 0 -0.4650 13.1180 -0.035 0.971723
## 1991 - 1989 == 0 15.3136 8.9695 1.707 0.087768 .
## 1992 - 1989 == 0 9.6627 6.9592 1.388 0.164996
## 1993 - 1989 == 0 7.3151 7.3613 0.994 0.320356
## 1994 - 1989 == 0 13.0085 7.4668 1.742 0.081477 .
## 1995 - 1989 == 0 30.8916 7.9051 3.908 9.31e-05 ***
## 1996 - 1989 == 0 13.6036 7.8966 1.723 0.084940 .
## 1997 - 1989 == 0 8.3921 7.3180 1.147 0.251475
## 1998 - 1989 == 0 28.2541 13.3394 2.118 0.034167 *
## 1991 - 1990 == 0 15.7786 14.4506 1.092 0.274878
## 1992 - 1990 == 0 10.1277 13.1209 0.772 0.440190
## 1993 - 1990 == 0 7.7801 13.3694 0.582 0.560611
## 1994 - 1990 == 0 13.4735 13.3364 1.010 0.312361
## 1995 - 1990 == 0 31.3566 13.5850 2.308 0.020989 *
## 1996 - 1990 == 0 14.0686 13.7748 1.021 0.307099
## 1997 - 1990 == 0 8.8571 13.2657 0.668 0.504346
## 1998 - 1990 == 0 28.7191 17.5475 1.637 0.101704
## 1992 - 1991 == 0 -5.6509 9.0140 -0.627 0.530724
## 1993 - 1991 == 0 -7.9984 9.1313 -0.876 0.381063
## 1994 - 1991 == 0 -2.3051 9.3109 -0.248 0.804471
## 1995 - 1991 == 0 15.5781 9.7907 1.591 0.111584
## 1996 - 1991 == 0 -1.7100 9.8262 -0.174 0.861846
## 1997 - 1991 == 0 -6.9215 8.6034 -0.805 0.421105
## 1998 - 1991 == 0 12.9405 14.7017 0.880 0.378748
## 1993 - 1992 == 0 -2.3475 7.5446 -0.311 0.755684
## 1994 - 1992 == 0 3.3459 7.4915 0.447 0.655152
## 1995 - 1992 == 0 21.2290 7.9271 2.678 0.007406 **
## 1996 - 1992 == 0 3.9409 7.9137 0.498 0.618495
## 1997 - 1992 == 0 -1.2706 7.3109 -0.174 0.862027
## 1998 - 1992 == 0 18.5914 13.3481 1.393 0.163675
## 1994 - 1993 == 0 5.6934 7.7298 0.737 0.461395
## 1995 - 1993 == 0 23.5765 8.4420 2.793 0.005226 **
## 1996 - 1993 == 0 6.2884 8.5219 0.738 0.460569
## 1997 - 1993 == 0 1.0769 7.6974 0.140 0.888732
## 1998 - 1993 == 0 20.9389 13.8495 1.512 0.130563
## 1995 - 1994 == 0 17.8831 8.1687 2.189 0.028581 *
## 1996 - 1994 == 0 0.5950 8.1562 0.073 0.941841
## 1997 - 1994 == 0 -4.6164 7.6059 -0.607 0.543880
## 1998 - 1994 == 0 15.2455 13.7648 1.108 0.268046
## 1996 - 1995 == 0 -17.2881 8.4289 -2.051 0.040261 *
## 1997 - 1995 == 0 -22.4996 8.0591 -2.792 0.005241 **
## 1998 - 1995 == 0 -2.6376 13.6383 -0.193 0.846649
## 1997 - 1996 == 0 -5.2115 8.1358 -0.641 0.521805
## 1998 - 1996 == 0 14.6505 13.6542 1.073 0.283287
## 1998 - 1997 == 0 19.8620 13.6881 1.451 0.146769
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## (Univariate p values reported)
Figure2<-grid.arrange(ENSO, Model1_plot, arrangeGrob(
Model2_plot, Model3b_plot, Model4b_plot, widths = c(3.9/12, 3.80/12, 4.30/12), ncol = 3),
nrow=3)
#ggsave(file="Fig1_V7.svg", plot=Figure1, width=6.69, height=6.69)
FIGURE 2 Illustration of the influence of most intense ENSO events on the coral cover in the ETP. (a) The upper panel indicates the positive (red) and negative (blue) SST anomalies in the Niño-3 region during 1971–2014 (NOAA Climate Prediction Center’s Extended Reconstructed Sea Surface Temperature ERSSTv5). The vertical orange bars indicate the strongest El Niño events. (b) The panel displays, at the horizontal axis, the 44-year live-coral cover aggregated at the region scale (n = 202) for ETP coral reefs with its fitted linear trend (gray dashed band). The fitted smooth line trend represents the long-term cycles of loss and recovery (blue line with the best polynomial fit with a 95% confidence level interval and a span = 0.2). The vertical axis shows the live coral cover. Box plots describe the minimum, maximum, interquartile ranges and median values of live-coral cover for each year. Notice the abrupt decrease in coral cover after 1982–83 El Niño; a smaller decline occurred after 1997–98. Random intercepts and random slopes models for the time-intervals: (c) 1974–1982, (d)1983–1997, and (e)1998–2014.
# Load packages
library(tidyverse)
library(dplyr)
library(data.table)
library(RColorBrewer)
# Select DHW data file: maxDHW.csv
dhw <- read.csv("Data/maxDHW.csv", header=TRUE)
# Exclude years 1981, and 2017-2019
dhw <- dhw[!dhw$Year == 1981, ]
# Create groups for plot
ts.Cano_island <- subset(dhw, dhw$Site == "Cano_island")
ts.Uva_island <- subset(dhw, dhw$Site == "Uva_Island")
ts.Gorgona_island <- subset(dhw, dhw$Site == "Gorgona_Island")
ts.Darwin_island <- subset(dhw, dhw$Site == "Darwin_North_Anchorage")
ts.Bahia_Banderas <- subset(dhw, dhw$Site == "Bahia_Banderas")
# Year, maxDHW
par(oma = c(1,1,1,1))
par(mar = c(4,4,1,1), cex.axis = 1, las = 1)
par(lwd = 1)
plot(ts.Cano_island$Year,ts.Cano_island$maxDHW, # Create a white false line with real data
xlab="Years", ylab="Degree heating weeks (°C-weeks)", cex=1.5, cex.lab=0.6, col="white",axes = FALSE,
ylim = c(0, 28), lwd = 1.5)
par(new=T)
rect(1982, 0, 1983.5, 28, col="#fee8c8", border ="#fee8c8") # Strong ENSO years
par(new=T)
rect(1997, 0, 1998.5, 28, col="#fee8c8", border ="#fee8c8")
par(new=T)
plot(ts.Cano_island$Year, ts.Cano_island$maxDHW,type="l", # Thermally stable
xlab="", ylab="", lwd = 1.5, col="#8c2d04",axes = FALSE,
ylim = c(0, 28))
par(new=T)
plot(ts.Uva_island$Year, ts.Uva_island$maxDHW, type="l", # Thermally stable
xlab="", ylab="", lwd = 1.5, col="#feb24c",axes = FALSE,
ylim = c(0, 28))
par(new=T)
plot(ts.Gorgona_island$Year,ts.Gorgona_island$maxDHW, type="l", # Tropical upwelling
xlab="", ylab="",lwd = 1.5, col="#fb6a4a", axes = FALSE,
ylim = c(0, 28))
par(new=T)
plot(ts.Darwin_island$Year,ts.Darwin_island$maxDHW, type="l", # Galapagos
xlab="", ylab="", lwd = 1.5, col="#9e9ac8",axes = FALSE,
ylim = c(0, 28))
par(new=T)
plot(ts.Bahia_Banderas$Year, ts.Bahia_Banderas$maxDHW, type="l", # Seasonal
xlab="", ylab="",lwd = 1.5,col="#a6d96a", axes = FALSE,
ylim = c(0, 28))
y <- c(0,4,8,12,16,20,24,28) # Axis
axis(2,lwd = 0.5, at = y,cex.axis=0.6)
box(lwd = 0.5)
x <-c(1982:2014)
YearsR = as.character(c(1982:2014))
axis(1, at = x, labels = YearsR, lwd = 0.4, cex.axis=0.5, las = 2)
abline(h=4, col="#525252", lty=2)
abline(h=8, col="#525252", lty=2)
# Legend
sites <- c("Thermally stable - Caño Island","Thermally stable - Uva Island",
"Tropical upwelling - Gorgona Island","Equatorial Upwelling - Darwin Island",
"Seasonal - Bahía Banderas")
par(xpd=TRUE)
legend(x = 2000, y = 20,
xpd = NA, # or TRUE
inset=c(-0.5,0),
cex = 0.4,
legend = sites, # box categories
col = c("#8c2d04","#feb24c","#fb6a4a","#9e9ac8","#a6d96a"),
lty= 1,
lwd = 1.5,
x.intersp = 0.5,
y.intersp = 2,
bty="n")
We tested the hypothesis that the strength of the Log_annual_rate_of_change (response variable) is a function of the maxDHW experienced. We had measured the rate of change from four Thermal regimes, and fitted Thermal regimes as a random intercept, and estimate a common slope (change in the rate of change) for maxDHW across all Thermal regimes by fitting it as a fixed effect.
## Load data
roc.cover <- read.csv("Data/RoC_DHW.csv", header=TRUE) ### Select file: RoC_DWH.csv
# Load packages
library(lme4)
library(lmerTest)
library(MuMIn)
model6 <- lmerTest::lmer(Log_annual_rate_of_change ~ maxDHW + (1|Thermal_regime), data=roc.cover)
step(model6)
## Backward reduced random-effect table:
##
## Eliminated npar logLik AIC LRT Df Pr(>Chisq)
## <none> 4 -17.738 43.476
## (1 | Thermal_regime) 0 3 -27.772 61.545 20.069 1 7.469e-06
##
## <none>
## (1 | Thermal_regime) ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Backward reduced fixed-effect table:
## Degrees of freedom method: Satterthwaite
##
## Eliminated Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## maxDHW 0 0.65054 0.65054 1 128.59 9.873 0.002082 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Model found:
## Log_annual_rate_of_change ~ maxDHW + (1 | Thermal_regime)
summary(model6)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Log_annual_rate_of_change ~ maxDHW + (1 | Thermal_regime)
## Data: roc.cover
##
## REML criterion at convergence: 35.5
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -4.7962 -0.1472 0.1167 0.4152 1.8621
##
## Random effects:
## Groups Name Variance Std.Dev.
## Thermal_regime (Intercept) 0.02501 0.1581
## Residual 0.06589 0.2567
## Number of obs: 133, groups: Thermal_regime, 4
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) -0.085763 0.084256 3.155381 -1.018 0.38032
## maxDHW -0.016316 0.005193 128.594887 -3.142 0.00208 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr)
## maxDHW -0.165
r.squaredGLMM(model6) # returns the marginal and the conditional R2
## R2m R2c
## [1,] 0.05249512 0.3131505
# marginal R2 describes the proportion of variance explained by the fixed factor(s) alone:
# conditional R2 describes the proportion of variance explained by both the fixed and random factors
# Nakagawa, S. & Schielzeth, H. (2013).
# Load packages
library(tidyverse)
library(dplyr)
library(data.table)
library(RColorBrewer)
# Create groups for plot
Seasonal = subset(roc.cover, roc.cover$Thermal_regime == "Seasonal")
Thermally_Stable = subset(roc.cover, roc.cover$Thermal_regime == "Thermally_Stable")
Tropical_Upwelling = subset(roc.cover, roc.cover$Thermal_regime == "Upwelling")
Galapagos = subset(roc.cover, roc.cover$Thermal_regime == "Galapagos_Islands")
# Plot
par(oma = c(1,1,1,1))
par(mar = c(4,4,1,1), cex.axis = 1, las = 1)
par(lwd = 1)
plot (Seasonal$maxDHW, Seasonal$Log_annual_rate_of_change,
pch = 21, lwd = 0.5, bg = "#a6d96a93", col="black", axes = FALSE,
cex=1.5,cex.lab=0.6,
xlim = c(0,30), ylim = c(-1.75, 1),
xlab = "Degree heating weeks (°C-weeks)",
ylab = "Annual rate of change in coral cover")
par(new=T);
plot (Thermally_Stable$maxDHW, Thermally_Stable$Log_annual_rate_of_change,
xlab="", ylab="", pch = 21, lwd = 0.5, col="black", bg = "#feb24c93", axes = FALSE,
cex=1.5,cex.lab=0.6,
xlim = c(0,30), ylim = c(-1.75, 1))
par(new=T);
plot (Tropical_Upwelling$maxDHW, Tropical_Upwelling$Log_annual_rate_of_change,
xlab="", ylab="", pch = 21, lwd = 0.5, bg = "#fb6a4a93",col="black", axes = FALSE,
cex=1.5,cex.lab=0.6,
xlim = c(0,30), ylim = c(-1.75, 1))
par(new=T);
plot (Galapagos$maxDHW, Galapagos$Log_annual_rate_of_change,
xlab="", ylab="", pch = 21, lwd = 0.5, bg = "#9e9ac893",col="black", axes = FALSE,
cex=1.5,cex.lab=0.6,
xlim = c(0,30), ylim = c(-1.75, 1))
box(lwd = 0.5)
x <-c(0,5,10,15,20,25,30)
axis(1, lwd = 0.5, at = x, cex.axis = 0.5)
y <- c(-1.5,-1.0,-0.5,0,0.5,1)
axis(2,lwd = 0.5, at = y, cex.axis = 0.5)
abline(v=4, col="#525252", lty=2); abline(v=8, col="#525252", lty=2); abline(h=0, col="#525252", lty=2);
# Axis text I, II, III, IV
text(3.5, 0.3, "II", cex = 0.6); text(5, 0.3, "I", cex = 0.6)
text(3.5, - 0.3, "III", cex = 0.6); text(5, -0.3, "IV", cex = 0.6)
par(xpd=TRUE)
legend(x = 10, y = -0.5,
xpd = NA,
inset=c(-0.5,0),
legend = c("Equatorial Upwelling", "Thermally stable", "Tropical upwelling","Seasonal"),
pch=21,
pt.bg = c("#9e9ac8","#fb6a4a","#fb6a4a","#a6d96a"),
x.intersp = 0.5,
y.intersp = 2,
bty="n")
FIGURE 3 Heat stress and its effect on coral cover across the ETP. (a) The sequence of yearly maximum DHW experienced by the best-represented individual coral reefs sites from 1982 to 2014. When cumulative thermal stress reach 4°C-weeks and 8°C-weeks or higher (horizontal dotted lines), bleaching is likely, as well as coral mortality from thermal stress (Liu et al., 2006), respectively. Although the thresholds were developed in seasonal systems, and are not universally applicable for equatorial regions, they still useful for separating responses into groups. (b) The association between the maximum DHW and the coral cover change at each of four thermal regimes (dot colors as Fig. 2). Each point (n = 133) represents the DHW and the coral cover annual rate of change per individual coral reef site. The horizontal dotted line divides the positive (> 0) and the negative (< 0) rate of change caused by coral growth or mortality, respectively. The vertical dotted line indicates the DHW thresholds of 4°C-weeks and 8°C-weeks, respectively. The vertical and horizontal dotted lines create four response quadrants: I) heat stress and positive coral cover change (resilient response), II) no heat stress and positive coral cover change, III) no heat stress and negative coral cover change, and IV) heat stress and negative coral cover change.
No stress
no.stress <- subset(roc.cover, maxDHW < 4) # subset of DHW < 4
length(no.stress$maxDHW)/length(roc.cover$maxDHW) # ratio no stress: total
## [1] 0.8496241
1- length(no.stress$maxDHW)/length(roc.cover$maxDHW) # ratio stress: total
## [1] 0.1503759
positive.rate <- subset(no.stress, Log_annual_rate_of_change > 0) # subset of Rate > 0
length(positive.rate$Log_annual_rate_of_change) # Number of observations in
## [1] 48
negative.rate <- subset(no.stress, Log_annual_rate_of_change < 0) #subset of Rate > 0
length(negative.rate$Log_annual_rate_of_change) # Number of observations in
## [1] 58
No stress seasonal
no.stress.seasonal <- subset(no.stress, Thermal_regime == "Seasonal")
positive.rate <- subset(no.stress.seasonal, Log_annual_rate_of_change > 0)
negative.rate <- subset(no.stress.seasonal, Log_annual_rate_of_change < 0)
length(positive.rate$Log_annual_rate_of_change)/length(negative.rate$Log_annual_rate_of_change)
## [1] 0.07142857
No stress upwelling
no.stress.upwelling <- subset(no.stress, Thermal_regime == "Upwelling")
positive.rate <- subset(no.stress.upwelling, Log_annual_rate_of_change > 0)
negative.rate <- subset(no.stress.upwelling, Log_annual_rate_of_change < 0)
length(positive.rate$Log_annual_rate_of_change)/length(negative.rate$Log_annual_rate_of_change)
## [1] 1
No stress thermally stable
no.stress.Thermally_Stable <- subset(no.stress, Thermal_regime == "Thermally_Stable")
positive.rate <- subset(no.stress.Thermally_Stable, Log_annual_rate_of_change > 0)
negative.rate <- subset(no.stress.Thermally_Stable, Log_annual_rate_of_change < 0)
length(positive.rate$Log_annual_rate_of_change)/length(negative.rate$Log_annual_rate_of_change)
## [1] 1.2
No stress Galapagos
no.stress.Galapagos_Islands <- subset(no.stress, Thermal_regime == "Galapagos_Islands")
positive.rate <- subset(no.stress.Galapagos_Islands, Log_annual_rate_of_change > 0)
negative.rate <- subset(no.stress.Galapagos_Islands, Log_annual_rate_of_change < 0)
length(positive.rate$Log_annual_rate_of_change)/length(negative.rate$Log_annual_rate_of_change)
## [1] 1
library(pals)
library(dplyr)
library(reshape2)
library(RColorBrewer)
# Barplot number of surveys per year
obs_num = as.data.frame(table(cover$Year,cover$Region))
colnames(obs_num) = c("Year", "Region", "Freq")
wide_obs_num = dcast(obs_num, Region~Year, value = "Freq")
df = data.matrix(wide_obs_num[,2:45])
as.table(df)
## 1970 1971 1972 1974 1975 1976 1977 1978 1979 1980 1981 1982 1983 1984
## A 0 0 0 0 0 0 0 0 0 0 0 0 0 0
## B 0 0 0 0 0 0 0 0 0 0 0 0 0 0
## C 0 0 1 0 0 0 0 0 1 1 0 0 1 0
## D 0 0 0 0 0 0 0 0 0 0 0 0 0 0
## E 0 0 0 0 0 0 0 0 0 1 0 1 0 4
## F 0 0 0 2 3 2 0 0 0 0 0 1 0 1
## G 0 0 0 0 0 0 0 0 0 0 0 0 0 0
## H 3 10 2 13 4 4 3 5 2 5 2 3 5 3
## I 0 0 0 0 0 0 0 0 0 0 0 0 0 0
## J 0 0 0 0 0 0 0 0 0 0 0 0 0 0
## K 0 0 0 0 0 0 0 0 0 0 0 0 0 0
## L 0 0 0 0 0 0 0 2 0 0 0 0 0 0
## M 0 0 0 0 0 0 0 0 0 0 0 1 0 0
## N 2 6 0 2 1 0 0 0 0 0 0 1 4 2
## O 0 0 0 0 0 0 0 0 0 0 0 0 0 0
## P 0 0 0 0 0 0 0 0 0 0 0 0 0 0
## 1985 1986 1987 1988 1989 1990 1991 1992 1993 1994 1995 1996 1997 1998
## A 0 0 0 0 0 0 0 0 0 1 0 0 0 0
## B 0 0 5 0 0 0 0 1 0 0 0 0 0 0
## C 0 0 1 2 4 0 0 1 0 0 2 2 0 4
## D 0 0 0 0 0 0 0 0 0 0 0 0 0 0
## E 12 1 2 0 1 0 0 11 0 2 0 6 0 1
## F 0 1 1 0 1 0 0 0 2 0 0 0 0 0
## G 0 0 0 0 0 0 0 0 0 0 0 0 0 0
## H 4 5 2 3 1 1 0 4 1 4 4 0 1 1
## I 0 0 0 0 0 0 2 0 0 0 0 0 16 10
## J 0 0 1 1 1 0 3 1 1 0 0 0 2 2
## K 0 0 0 0 0 0 0 0 0 0 0 0 2 1
## L 0 0 0 0 0 0 0 0 0 1 1 4 4 2
## M 0 0 0 0 0 0 0 1 0 0 0 0 0 0
## N 0 1 5 1 2 0 1 1 1 1 1 1 1 1
## O 0 0 0 0 0 0 0 0 1 1 0 0 0 0
## P 0 0 0 0 0 0 0 0 0 0 0 0 0 0
## 1999 2000 2001 2002 2003 2004 2005 2006 2007 2008 2009 2010 2011 2012
## A 0 0 0 0 0 0 0 0 0 0 0 0 0 0
## B 0 0 0 1 0 0 0 1 0 0 4 0 0 0
## C 7 4 3 4 6 4 3 6 3 2 3 2 2 3
## D 0 0 0 0 0 0 0 0 0 0 0 4 0 1
## E 1 0 2 0 4 2 2 0 1 0 2 0 0 0
## F 0 0 1 1 0 1 3 6 2 2 0 0 1 0
## G 0 0 0 0 0 0 0 0 0 0 1 0 0 0
## H 0 1 4 14 12 1 1 1 1 0 1 2 0 0
## I 0 0 0 6 3 0 1 1 4 0 17 6 2 0
## J 0 0 0 1 2 0 0 0 0 5 7 0 0 1
## K 1 5 5 3 2 0 0 0 0 0 2 0 0 0
## L 0 0 2 0 1 0 1 3 2 0 17 3 1 1
## M 0 0 0 0 0 1 2 14 8 0 0 0 0 1
## N 1 1 1 1 1 1 0 0 0 0 1 0 0 0
## O 0 0 0 0 0 0 0 0 0 0 0 0 0 0
## P 0 0 0 0 0 2 0 0 0 0 0 0 0 0
## 2013 2014
## A 0 0
## B 1 1
## C 2 1
## D 4 2
## E 0 0
## F 0 0
## G 0 0
## H 0 0
## I 2 2
## J 0 0
## K 0 0
## L 1 3
## M 0 0
## N 0 0
## O 0 0
## P 0 0
df[ df>1 ] <- 1
levels(cover$Region) # How many regions? = how many colours
## [1] "Clipperton" "Cocos_Islands"
## [3] "Colombia" "Costa_Rica_Central"
## [5] "Costa_Rica_Southern" "Eastern_Galapagos_Islands"
## [7] "El_Salvador" "Gulf_of_Chiriqui"
## [9] "Mexican_Central" "Mexican_Northern"
## [11] "Mexican_Southern" "Nicaragua_Papagayo_Zone"
## [13] "Northern_Galapagos_Islands" "Panama_Bight"
## [15] "Revillagigedo_Islands" "Western_Galapagos_Islands"
Regions <- c(levels(obs_num$Region))
color.regions = warmcool(17)
# Figure 2a
par(mfcol=c(2,1), cex.main = 1.5, mar = c(4,4,1,1), cex.lab = 1.5, cex.axis = 1, las = 1)
barplot(df,
ylim=c(0,15),
col = color.regions,
ylab="Number of regions surveyed per year",
xlab="Year",
#lwd = 1,
las = 2,
axis.lty = 1,
cex.names = 0.4,
cex.axis = 0.8)
# Legend as a box
legend( #"bottomright",
x = 2, y = 16,
xpd = NA, # or TRUE
legend = Regions,
fill = color.regions[1:17],
inset=c(-0.5,0),
cex = 0.35,
text.font=0.5,
x.intersp = -0.5,
y.intersp = 1,
lwd = 0.1,
bty="n")
# Figure 2b
obs_num = as.data.frame(table(cover$Year,cover$Country)) # Table per country
colnames(obs_num) = c("Year", "Country", "Freq")
wide_obs_num = dcast(obs_num, Country~Year, value = "Freq")
dfc = data.matrix(wide_obs_num[,2:45])
as.table(dfc)
## 1970 1971 1972 1974 1975 1976 1977 1978 1979 1980 1981 1982 1983 1984
## A 0 0 1 0 0 0 0 0 1 1 0 0 1 0
## B 0 0 2 9 0 0 0 2 0 2 0 1 0 4
## C 0 0 0 2 3 2 0 0 0 0 0 2 0 1
## D 0 0 0 0 0 0 0 0 0 0 0 0 0 0
## E 0 0 0 0 0 0 0 0 0 0 0 0 0 0
## F 0 0 0 0 0 0 0 0 0 0 0 0 0 0
## G 0 0 0 0 0 0 0 0 0 0 0 0 0 0
## H 5 16 0 6 5 4 3 5 2 4 2 4 9 5
## 1985 1986 1987 1988 1989 1990 1991 1992 1993 1994 1995 1996 1997 1998
## A 0 0 1 2 4 0 0 1 0 0 2 2 0 4
## B 12 1 7 0 1 0 0 15 0 6 4 10 4 3
## C 0 1 1 0 1 0 0 1 2 0 0 0 0 0
## D 0 0 0 0 0 0 0 0 0 0 0 0 0 0
## E 0 0 0 0 0 0 0 0 0 1 0 0 0 0
## F 0 0 1 1 1 0 5 1 2 1 0 0 20 13
## G 0 0 0 0 0 0 0 0 0 0 0 0 0 0
## H 4 6 7 4 3 1 1 2 2 2 2 1 2 2
## 1999 2000 2001 2002 2003 2004 2005 2006 2007 2008 2009 2010 2011 2012
## A 7 4 3 4 6 4 3 6 3 2 3 2 2 3
## B 1 0 7 1 5 2 3 4 3 0 12 7 1 2
## C 0 0 1 1 0 4 5 20 10 2 0 0 1 1
## D 0 0 0 0 0 0 0 0 0 0 1 0 0 0
## E 0 0 0 0 0 0 0 0 0 0 0 0 0 0
## F 1 5 5 10 7 0 1 1 4 5 26 6 2 1
## G 0 0 0 0 0 0 0 0 0 0 11 0 0 0
## H 1 2 2 15 13 2 1 1 1 0 2 2 0 0
## 2013 2014
## A 2 1
## B 6 6
## C 0 0
## D 0 0
## E 0 0
## F 2 2
## G 0 0
## H 0 0
dfc[ dfc>1 ] <- 1 # replace all values greater than one with one
color.countries <- brewer.pal(9, "Oranges")
countries = c(levels(obs_num$Country))
barplot(dfc,
col = color.countries,
ylim=c(0,10),
ylab="Number of countries surveyed per year",
xlab="Year",
las = 2,
axis.lty = 1,
cex.names = 0.4,
cex.axis = 0.8)
# Legend as a box
legend( #"bottomright",
x = 2, y = 10,
xpd = NA, # or TRUE
legend = countries,
fill = color.countries[1:9],
inset=c(-0.5,0),
cex = 0.35,
text.font=1,
x.intersp = -0.5,
y.intersp = 1,
lwd = 0.1,
bty="n")
FIGURE S1 Regions and countries surveyed per year. (a) Number of regions surveyed per year in the ETP. Since 1970, the number of surveys in all regions increased by year up to a maximum of 10 regions surveyed in 2009. (b) Number of countries surveyed per year (Colombia, Costa Rica, Ecuador, El Salvador, France, Mexico, Nicaragua and Panama). Since 1970, the number of surveys in all the countries reached a maximum of 6 countries in 2009.
Temporal variation in the live coral cover for 1970–2014.
The smooth line represents the best fit to the time series with a 95% confidence level interval and a span = 0.3.
library(ggplot2)
library(lemon)
aggr.region_S2<-aggr.region
aggr.region_S2$Thermal_regime<-factor(aggr.region_S2$Thermal_regime,
levels=c("Tropical upwelling", "Seasonal",
"Thermally stable", "Equatorial upwelling"),
labels = c ("a", "b", "c", "d"))
FigureS2<-ggplot(aggr.region_S2, aes(x=Year, y=Coral_cover)) +
theme_bw() + theme(panel.border = element_blank(), # set ggplot to white background
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
axis.line = element_line(),
strip.background = element_blank(),
strip.text = element_text(face="bold",hjust = 0, vjust = -1))+
geom_point()+ geom_smooth (span = 0.3) +
scale_x_continuous(limits=c(1970, 2014),
breaks = seq(1970, 2014, by=5),
expand = c(0.01, 0.01), "Year") +
scale_y_continuous(limits=c(-19, 91), expand = c(0, 0), "Live coral cover (%)") +
facet_wrap(Thermal_regime~., ncol=1, scales = "free_x", strip.position ="top")
FigureS2
#ggsave(file="FigureS2.svg", plot=FigureS2, width=6.69, height=6.69)
FIGURE S2 Temporal variation in the live coral cover for 1970–2014. (a) Tropical upwelling regime (n = 75). (b) Seasonal regime (n = 28). (c) Thermally stable regime (n = 74). (d) Equatorial upwelling regime (n = 25). The smooth line represents the best fit to the time series with a 95% confidence level interval and a span = 0.3.
FigureS2b<-ggplot(cover, aes(x=Year, y=Coral_cover)) +
theme_bw() + theme(panel.border = element_blank(), # set ggplot to white background
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
axis.line = element_line(),
strip.background = element_blank(),
strip.text = element_text(face="bold",hjust = 0, vjust = -1))+
geom_point(aes(colour=Region)) +
geom_smooth (span = 0.3, colour="gray") +
scale_x_continuous(limits=c(1970, 2014),
breaks = seq(1970, 2014, by=5),
expand = c(0.01, 0.01), "Year") +
scale_y_continuous(limits=c(-19, 91), expand = c(0, 0), "Live coral cover (%)") +
facet_wrap(Thermal_regime~., ncol=1, scales = "free_x", strip.position ="top") +
theme(legend.position = "bottom", legend.direction = "horizontal" )
FigureS2b
We found in published studies four comprehensive monitoring datasets for Uva Island, Gorgona Island, Costa Rica, and the Eastern and Northern Galapagos Islands, which correspond to well recognized long-term monitoring programs
# Create groups
Uva_Island <- subset(aggr.site, Site == "Uva_Island")
Gorgona <-subset(aggr.site, Site == "La_Azufrada")
Cano_island <-subset(aggr.site, Site == "Cano_island")
EGI <-subset(aggr.site, Region == "Eastern_Galapagos_Islands")
NGI <-subset(aggr.site, Region == "Northern_Galapagos_Islands")
Monitoring_sites <- rbind(Uva_Island,Gorgona,Cano_island,EGI,NGI)
model7 <- lme(Coral_cover ~ -1 + Thermal_regime * Year, random = ~1|Region, data=Monitoring_sites)
summary(model7)
## Linear mixed-effects model fit by REML
## Data: Monitoring_sites
## AIC BIC logLik
## 978.1405 999.7444 -481.0703
##
## Random effects:
## Formula: ~1 | Region
## (Intercept) Residual
## StdDev: 14.54312 15.90812
##
## Fixed effects: Coral_cover ~ -1 + Thermal_regime * Year
## Value Std.Error DF t-value
## Thermal_regimeThermally stable -565.5374 449.5031 2 -1.2581391
## Thermal_regimeTropical upwelling 262.8581 757.1479 2 0.3471688
## Thermal_regimeEquatorial upwelling 373.7671 418.4179 2 0.8932865
## Year 0.2928 0.2257 109 1.2971984
## Thermal_regimeTropical upwelling:Year -0.3987 0.4404 109 -0.9054869
## Thermal_regimeEquatorial upwelling:Year -0.4645 0.3078 109 -1.5092165
## p-value
## Thermal_regimeThermally stable 0.3353
## Thermal_regimeTropical upwelling 0.7616
## Thermal_regimeEquatorial upwelling 0.4660
## Year 0.1973
## Thermal_regimeTropical upwelling:Year 0.3672
## Thermal_regimeEquatorial upwelling:Year 0.1341
## Correlation:
## Thr_Ts Thr_Tu Thr_Eu Year T_Tu:Y
## Thermal_regimeTropical upwelling 0.000
## Thermal_regimeEquatorial upwelling 0.000 0.000
## Year -1.000 0.000 0.000
## Thermal_regimeTropical upwelling:Year 0.512 -0.858 0.000 -0.513
## Thermal_regimeEquatorial upwelling:Year 0.733 0.000 -0.680 -0.733 0.376
##
## Standardized Within-Group Residuals:
## Min Q1 Med Q3 Max
## -2.73438942 -0.56940133 -0.09515832 0.50988895 3.13690545
##
## Number of Observations: 116
## Number of Groups: 5
anova(model7)
## numDF denDF F-value p-value
## Thermal_regime 3 2 7.535454 0.1194
## Year 1 109 0.023766 0.8778
## Thermal_regime:Year 2 109 1.205466 0.3035
plot(ranef(model7))
plot(model7)
resnorm7 <- resid(model7)
hist(resnorm7, xlab = "Residuals", main = "")
coef.m7 <- as.data.frame(coef(summary(model7)))
plot(model7, resid(., type = "p") ~ fitted(.) | Region, abline = 0)
plot(model7, Coral_cover ~ fitted(.) | Region, abline = c(0,1))
plot(model7, Coral_cover ~ fitted(.) | Thermal_regime, abline = c(0,1))
FigureS3 <-ggplot(Monitoring_sites, aes(x=Year, y=Coral_cover)) +
geom_boxplot(aes(group=Year), outlier.shape = NA) +
stat_boxplot(aes(group=Year), geom = 'errorbar')+
geom_smooth(method = lm, se=FALSE, linetype = "dashed", colour="red")+
geom_smooth(span = 0.2, se=T, colour="darkgray")+
scale_y_continuous("Live coral cover (%)", limits = c(-15, 95),
breaks = seq(0, 90, by=20), expand = c(0,0))+
scale_x_continuous("", limits = c(1969, 2015),
breaks = seq(1970, 2014, by=2), expand = c(0,0))+
annotate("rect", xmin = 1982, xmax = 1983, ymin = 0, ymax = 80, alpha = .2, fill="red")+
annotate("rect", xmin = 1997, xmax = 1998, ymin = 0, ymax = 80, alpha = .2, fill="red")+
annotate("rect", xmin = 1972, xmax = 1973, ymin = 0, ymax = 80, alpha = .2, fill="red")
FigureS3
FigureS3b<- FigureS3 + geom_point(aes(colour=Location), alpha=0.4) +
theme(legend.position = "bottom")
FigureS3b
#ggsave(file="FigS3.png", plot=FigureS3, width=6.69, height=4)
#ggsave(file="FigS3b.png", plot=FigureS3b, width=6.69, height=4.69)
Model 1: All years (1970_2014) removing each termal regime
Model1_NoEquatorial_Upwelling <- ggplot(data=subset(aggr.region, Thermal_regime!="Tropical upwelling"), aes(x=Year, y=Coral_cover)) +
geom_boxplot(aes(group=Year), outlier.shape = NA) +
stat_boxplot(aes(group=Year), geom = 'errorbar')+
geom_smooth(method = lm, se=FALSE, linetype = "dashed", colour="red")+
geom_smooth(span = 0.2, se=T, colour="darkgray")+
#geom_point(aes(colour=Thermal_regime), alpha=0.5) +
scale_y_continuous("Live coral cover (%)", limits = c(-5, 90),
breaks = seq(0, 90, by=20), expand = c(0,0))+
scale_x_continuous("", limits = c(1969, 2015),
breaks = seq(1970, 2014, by=2), expand = c(0,0))+
scale_color_manual(values=my_colours) + labs(title="a")+
#annotate("text", x = 1983, y = 70, label = "*", size=8)+
annotate("rect", xmin = 1982, xmax = 1983, ymin = 0, ymax = 80, alpha = .2, fill="red")+
annotate("rect", xmin = 1997, xmax = 1998, ymin = 0, ymax = 80, alpha = .2, fill="red")+
annotate("rect", xmin = 1972, xmax = 1973, ymin = 0, ymax = 80, alpha = .2, fill="red")
#Model1_NoEquatorial_Upwelling
Model1_NoSeasonal <- ggplot(data=subset(aggr.region, Thermal_regime!="Seasonal"), aes(x=Year, y=Coral_cover)) +
geom_boxplot(aes(group=Year), outlier.shape = NA) +
stat_boxplot(aes(group=Year), geom = 'errorbar')+
geom_smooth(method = lm, se=FALSE, linetype = "dashed", colour="red")+
geom_smooth(span = 0.2, se=T, colour="darkgray")+
#geom_point(aes(colour=Thermal_regime), alpha=0.5) +
scale_y_continuous("Live coral cover (%)", limits = c(-5, 90),
breaks = seq(0, 90, by=20), expand = c(0,0))+
scale_x_continuous("", limits = c(1969, 2015),
breaks = seq(1970, 2014, by=2), expand = c(0,0))+
scale_color_manual(values=my_colours) + labs(title="b")+
#annotate("text", x = 1983, y = 70, label = "*", size=8)+
annotate("rect", xmin = 1982, xmax = 1983, ymin = 0, ymax = 80, alpha = .2, fill="red")+
annotate("rect", xmin = 1997, xmax = 1998, ymin = 0, ymax = 80, alpha = .2, fill="red")+
annotate("rect", xmin = 1972, xmax = 1973, ymin = 0, ymax = 80, alpha = .2, fill="red")
#Model1_NoSeasonal
Model1_NoThermallyStable <- ggplot(data=subset(aggr.region, Thermal_regime!="Thermally stable"), aes(x=Year, y=Coral_cover)) +
geom_boxplot(aes(group=Year), outlier.shape = NA) +
stat_boxplot(aes(group=Year), geom = 'errorbar')+
geom_smooth(method = lm, se=FALSE, linetype = "dashed", colour="red")+
geom_smooth(span = 0.2, se=T, colour="darkgray")+
#geom_point(aes(colour=Thermal_regime), alpha=0.5) +
scale_y_continuous("Live coral cover (%)", limits = c(-5, 90),
breaks = seq(0, 90, by=20), expand = c(0,0))+
scale_x_continuous("", limits = c(1969, 2015),
breaks = seq(1970, 2014, by=2), expand = c(0,0))+
scale_color_manual(values=my_colours) + labs(title="c")+
#annotate("text", x = 1983, y = 70, label = "*", size=8)+
annotate("rect", xmin = 1982, xmax = 1983, ymin = 0, ymax = 80, alpha = .2, fill="red")+
annotate("rect", xmin = 1997, xmax = 1998, ymin = 0, ymax = 80, alpha = .2, fill="red")+
annotate("rect", xmin = 1972, xmax = 1973, ymin = 0, ymax = 80, alpha = .2, fill="red")
# Model1_NoThermallyStable
Model1_No_EqUpwelling <- ggplot(data=subset(aggr.region, Thermal_regime!="Equatorial upwelling"), aes(x=Year, y=Coral_cover)) +
geom_boxplot(aes(group=Year), outlier.shape = NA) +
stat_boxplot(aes(group=Year), geom = 'errorbar')+
geom_smooth(method = lm, se=FALSE, linetype = "dashed", colour="red")+
geom_smooth(span = 0.2, se=T, colour="darkgray")+
#geom_point(aes(colour=Thermal_regime), alpha=0.5) +
scale_y_continuous("Live coral cover (%)", limits = c(-5, 90),
breaks = seq(0, 90, by=20), expand = c(0,0))+
scale_x_continuous("", limits = c(1969, 2015),
breaks = seq(1970, 2014, by=2), expand = c(0,0))+
scale_color_manual(values=my_colours) + labs(title="d")+
#annotate("text", x = 1983, y = 70, label = "*", size=8)+
annotate("rect", xmin = 1982, xmax = 1983, ymin = 0, ymax = 80, alpha = .2, fill="red")+
annotate("rect", xmin = 1997, xmax = 1998, ymin = 0, ymax = 80, alpha = .2, fill="red")+
annotate("rect", xmin = 1972, xmax = 1973, ymin = 0, ymax = 80, alpha = .2, fill="red")
# Model1_No_EqUpwelling
FigureS4<-grid.arrange(Model1_NoEquatorial_Upwelling, Model1_NoSeasonal,
Model1_NoThermallyStable, Model1_No_EqUpwelling,
ncol = 2)
#ggsave(file="FigS5.png", plot=FigureS5, width=6.69, height=6.69)
#ggsave(file="FigS5.pdf", plot=FigureS5, width=6.69, height=6.69)
Figure S4: Coral cover trends excluding (a) Thermally stable, (b) Tropical Upwelling, (c) Equatorial upwelling, and (d) Seasonal thermal regimes
To use “bootMer” models need to be run with lme4 instead nlme. Please check accuracy of the new (b) models
# 1. Build lme4 model:
model1b <- lmer(
Coral_cover ~ Thermal_regime * Year +
(1|Region), data=aggr.region, REML = FALSE)
summary(model1b)
## Linear mixed model fit by maximum likelihood . t-tests use
## Satterthwaite's method [lmerModLmerTest]
## Formula: Coral_cover ~ Thermal_regime * Year + (1 | Region)
## Data: aggr.region
##
## AIC BIC logLik deviance df.resid
## 1746.3 1779.4 -863.2 1726.3 192
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.32243 -0.72526 0.00017 0.57206 3.08905
##
## Random effects:
## Groups Name Variance Std.Dev.
## Region (Intercept) 96.1 9.803
## Residual 268.7 16.393
## Number of obs: 202, groups: Region, 16
##
## Fixed effects:
## Estimate Std. Error df
## (Intercept) -425.34274 371.24038 200.81720
## Thermal_regimeTropical upwelling 155.93891 520.89698 201.08927
## Thermal_regimeEquatorial upwelling 661.72680 677.59229 195.76846
## Thermal_regimeSeasonal 1436.94164 901.75966 201.77989
## Year 0.22505 0.18575 200.56729
## Thermal_regimeTropical upwelling:Year -0.07308 0.26052 200.89738
## Thermal_regimeEquatorial upwelling:Year -0.33411 0.33909 195.36095
## Thermal_regimeSeasonal:Year -0.71540 0.45058 201.68330
## t value Pr(>|t|)
## (Intercept) -1.146 0.253
## Thermal_regimeTropical upwelling 0.299 0.765
## Thermal_regimeEquatorial upwelling 0.977 0.330
## Thermal_regimeSeasonal 1.593 0.113
## Year 1.212 0.227
## Thermal_regimeTropical upwelling:Year -0.280 0.779
## Thermal_regimeEquatorial upwelling:Year -0.985 0.326
## Thermal_regimeSeasonal:Year -1.588 0.114
##
## Correlation of Fixed Effects:
## (Intr) Thr_Tu Thr_Eu Thrm_S Year T_Tu:Y T_Eu:Y
## Thrml_rgmTu -0.713
## Thrml_rgmEu -0.548 0.390
## Thrml_rgmSs -0.397 0.283 0.218
## Year -1.000 0.713 0.548 0.397
## Thrml_rTu:Y 0.713 -1.000 -0.391 -0.283 -0.713
## Thrml_rEu:Y 0.548 -0.390 -1.000 -0.217 -0.548 0.391
## Thrml_rgS:Y 0.398 -0.283 -0.218 -1.000 -0.397 0.283 0.217
anova(model1b) # Thermal regime is not significant anymore??
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Thermal_regime 839.67 279.889 3 199.65 1.0415 0.3753
## Year 41.27 41.273 1 199.94 0.1536 0.6956
## Thermal_regime:Year 843.71 281.237 3 199.47 1.0465 0.3731
ranova(model1b) # Region seems the only significant effect
## ANOVA-like table for random-effects: Single term deletions
##
## Model:
## Coral_cover ~ Thermal_regime + Year + (1 | Region) + Thermal_regime:Year
## npar logLik AIC LRT Df Pr(>Chisq)
## <none> 10 -863.15 1746.3
## (1 | Region) 9 -870.85 1759.7 15.399 1 8.705e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
drop1(model1b)
## Single term deletions using Satterthwaite's method:
##
## Model:
## Coral_cover ~ Thermal_regime * Year + (1 | Region)
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Thermal_regime:Year 843.71 281.24 3 199.47 1.0465 0.3731
Not sure if this if right since the model itself is not significant
# 2. Predict values:
pred1 <- predict(model1b,re.form = NA)
#3. Bootstrap CI:
boot1 <- bootMer(model1b, predict, nsim = 1000, re.form = NULL) # include random effects, reduce CI lot!
std.err <- apply(boot1$t, 2, sd)
CI.lo_1 <- pred1 - std.err*1.96
CI.hi_1 <- pred1 + std.err*1.96
#Plot
Model1b_plot<- ggplot(
aggr.region, aes(x = Year, y = Coral_cover, colour = Thermal_regime)) +
geom_line(aes(y = pred1),size=2) +
geom_point(aes(fill=factor(Thermal_regime)),
shape = 21, colour = "black", size = 2, stroke = 0.3, alpha=0.5) +
geom_ribbon(aes(ymin = CI.lo_1, ymax = CI.hi_1),
size=2, alpha = 0.03, linetype = 0) +
scale_color_manual(values=my_colours) +
scale_fill_manual(values=my_colours) +
scale_y_continuous("Live coral cover (%)",
limits = c(-20, 90),
breaks = seq(0, 90, by=20), expand = c(0,0))+
scale_x_continuous("", limits = c(1969, 2015),
breaks = seq(1970, 2014, by=2), expand = c(0,0))
Model1b_plot
Not enough data
# 1. Build lme4 model:
model3b <- lmer(
Coral_cover ~ Thermal_regime * Year +
(1|Region), data=ENSO.83_97, REML = FALSE)
summary(model3b)
## Linear mixed model fit by maximum likelihood . t-tests use
## Satterthwaite's method [lmerModLmerTest]
## Formula: Coral_cover ~ Thermal_regime * Year + (1 | Region)
## Data: ENSO.83_97
##
## AIC BIC logLik deviance df.resid
## 499.0 520.8 -239.5 479.0 55
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.0243 -0.5395 -0.0293 0.5949 3.2102
##
## Random effects:
## Groups Name Variance Std.Dev.
## Region (Intercept) 7.645 2.765
## Residual 86.995 9.327
## Number of obs: 65, groups: Region, 10
##
## Fixed effects:
## Estimate Std. Error df
## (Intercept) -2558.7536 885.2036 59.4204
## Thermal_regimeTropical upwelling -4663.3240 1215.3575 62.0348
## Thermal_regimeEquatorial upwelling 2670.5406 2851.0741 56.7512
## Thermal_regimeSeasonal 8827.9001 2025.6205 62.6391
## Year 1.2929 0.4449 59.3736
## Thermal_regimeTropical upwelling:Year 2.3555 0.6105 61.9948
## Thermal_regimeEquatorial upwelling:Year -1.3488 1.4341 56.7453
## Thermal_regimeSeasonal:Year -4.4231 1.0170 62.5994
## t value Pr(>|t|)
## (Intercept) -2.891 0.005362 **
## Thermal_regimeTropical upwelling -3.837 0.000294 ***
## Thermal_regimeEquatorial upwelling 0.937 0.352893
## Thermal_regimeSeasonal 4.358 4.97e-05 ***
## Year 2.906 0.005132 **
## Thermal_regimeTropical upwelling:Year 3.858 0.000275 ***
## Thermal_regimeEquatorial upwelling:Year -0.941 0.350946
## Thermal_regimeSeasonal:Year -4.349 5.13e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) Thr_Tu Thr_Eu Thrm_S Year T_Tu:Y T_Eu:Y
## Thrml_rgmTu -0.728
## Thrml_rgmEu -0.310 0.226
## Thrml_rgmSs -0.430 0.313 0.133
## Year -1.000 0.728 0.310 0.430
## Thrml_rTu:Y 0.729 -1.000 -0.226 -0.313 -0.729
## Thrml_rEu:Y 0.310 -0.226 -1.000 -0.133 -0.310 0.226
## Thrml_rgS:Y 0.430 -0.313 -0.133 -1.000 -0.430 0.313 0.133
anova(model3b) # Thermal regime and its interaciton with year are significant
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Thermal_regime 4445.7 1481.89 3 60.690 17.0343 3.858e-08 ***
## Year 87.0 86.98 1 60.112 0.9998 0.3214
## Thermal_regime:Year 4462.4 1487.47 3 60.665 17.0984 3.672e-08 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
ranova(model3b) # ramdom effects are not
## ANOVA-like table for random-effects: Single term deletions
##
## Model:
## Coral_cover ~ Thermal_regime + Year + (1 | Region) + Thermal_regime:Year
## npar logLik AIC LRT Df Pr(>Chisq)
## <none> 10 -239.52 499.03
## (1 | Region) 9 -240.10 498.20 1.1616 1 0.2811
drop1(model3b)
## Single term deletions using Satterthwaite's method:
##
## Model:
## Coral_cover ~ Thermal_regime * Year + (1 | Region)
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Thermal_regime:Year 4462.4 1487.5 3 60.665 17.098 3.672e-08 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# 2. Predict values:
pred3 <- predict(model3b,re.form = NA)
#3. Bootstrap CI:
boot3 <- bootMer(model3b, predict, nsim = 1000, re.form = NA) # No ramndom effects
std.err <- apply(boot3$t, 2, sd)
CI.lo_3 <- pred3 - std.err*1.96
CI.hi_3 <- pred3 + std.err*1.96
#Plot
Model3b_plot<- ggplot(
ENSO.83_97, aes(x = Year, y = Coral_cover, colour = Thermal_regime))+
geom_line(aes(y = pred3),size=2) +
geom_point(aes(fill=factor(Thermal_regime)),
shape = 21, colour = "black", size = 2, stroke = 0.3, alpha=0.5) +
geom_ribbon(aes(ymin = CI.lo_3, ymax = CI.hi_3),
size=2, alpha = 0.03, linetype = 0) +
scale_fill_manual(values=my_colours)+
scale_colour_manual(values=my_colours) +
scale_x_continuous("", limits = c(1982.7, 1998.2),
breaks = seq(1983, 1997, by=2),
expand = c(0.02,0.02))+
scale_y_continuous("Live coral cover (%)",
limits = c(-20, 90),
breaks = seq(0, 90, by=20), expand = c(0,0))
Model3b_plot
#ggsave(file="Model3b_plot.png", plot=Model3b_plot, width=4, height=4)
#ggsave(file="Model3b_plot.pdf", plot=Model3b_plot, width=4, height=4)
# 1. Build lme4 model:
model4b <- lmer(
Coral_cover ~ Thermal_regime * Year +
(1|Region), data=ENSO.98_14, REML = FALSE)
summary(model4b)
## Linear mixed model fit by maximum likelihood . t-tests use
## Satterthwaite's method [lmerModLmerTest]
## Formula: Coral_cover ~ Thermal_regime * Year + (1 | Region)
## Data: ENSO.98_14
##
## AIC BIC logLik deviance df.resid
## 894.7 921.2 -437.4 874.7 94
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.97345 -0.54963 -0.01987 0.30764 2.79401
##
## Random effects:
## Groups Name Variance Std.Dev.
## Region (Intercept) 81.52 9.029
## Residual 224.67 14.989
## Number of obs: 104, groups: Region, 12
##
## Fixed effects:
## Estimate Std. Error df
## (Intercept) -1734.4099 1338.8526 103.8926
## Thermal_regimeTropical upwelling 4885.3253 1681.7337 103.9991
## Thermal_regimeEquatorial upwelling -332.4797 3113.7676 96.2513
## Thermal_regimeSeasonal 1970.4798 2034.4368 100.1488
## Year 0.8782 0.6677 103.9042
## Thermal_regimeTropical upwelling:Year -2.4312 0.8385 103.9997
## Thermal_regimeEquatorial upwelling:Year 0.1641 1.5523 96.2303
## Thermal_regimeSeasonal:Year -0.9798 1.0140 100.0511
## t value Pr(>|t|)
## (Intercept) -1.295 0.19804
## Thermal_regimeTropical upwelling 2.905 0.00449 **
## Thermal_regimeEquatorial upwelling -0.107 0.91519
## Thermal_regimeSeasonal 0.969 0.33510
## Year 1.315 0.19132
## Thermal_regimeTropical upwelling:Year -2.899 0.00456 **
## Thermal_regimeEquatorial upwelling:Year 0.106 0.91605
## Thermal_regimeSeasonal:Year -0.966 0.33623
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) Thr_Tu Thr_Eu Thrm_S Year T_Tu:Y T_Eu:Y
## Thrml_rgmTu -0.796
## Thrml_rgmEu -0.430 0.342
## Thrml_rgmSs -0.652 0.519 0.280
## Year -1.000 0.796 0.430 0.652
## Thrml_rTu:Y 0.796 -1.000 -0.342 -0.519 -0.796
## Thrml_rEu:Y 0.430 -0.342 -1.000 -0.280 -0.430 0.343
## Thrml_rgS:Y 0.652 -0.519 -0.280 -1.000 -0.652 0.519 0.280
anova(model4b) # Thermal regime and its interaciton with year are significant
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Thermal_regime 2257.18 752.39 3 99.562 3.3489 0.02210 *
## Year 4.86 4.86 1 96.337 0.0217 0.88333
## Thermal_regime:Year 2248.88 749.63 3 99.519 3.3366 0.02244 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
ranova(model4b) # ramdom effects are significant
## ANOVA-like table for random-effects: Single term deletions
##
## Model:
## Coral_cover ~ Thermal_regime + Year + (1 | Region) + Thermal_regime:Year
## npar logLik AIC LRT Df Pr(>Chisq)
## <none> 10 -437.36 894.73
## (1 | Region) 9 -443.08 904.15 11.42 1 0.0007265 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
drop1(model4b)
## Single term deletions using Satterthwaite's method:
##
## Model:
## Coral_cover ~ Thermal_regime * Year + (1 | Region)
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Thermal_regime:Year 2248.9 749.63 3 99.519 3.3366 0.02244 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# 2. Predict values:
pred4 <- predict(model4b,re.form = NA)
#3. Bootstrap CI:
boot4 <- bootMer(model4b, predict, nsim = 1000, re.form = NULL) # Include ramndom effects
std.err <- apply(boot4$t, 2, sd)
CI.lo_4 <- pred4 - std.err*1.96
CI.hi_4 <- pred4 + std.err*1.96
#Plot
Model4b_plot<- ggplot(ENSO.98_14, aes(
x = Year, y = Coral_cover, colour = Thermal_regime)) +
geom_line(aes(y = pred4),size=2) +
geom_point(aes(fill=factor(Thermal_regime)),
shape = 21, colour = "black", size = 2,
stroke = 0.3, alpha=0.5) +
geom_ribbon(aes(ymin = CI.lo_4, ymax = CI.hi_4),size=2,
alpha = 0.03, linetype = 0) +
scale_y_continuous("Live coral cover (%)",
limits = c(-20, 90),
breaks = seq(0, 90, by=20), expand = c(0,0))+
scale_x_continuous("", limits = c(1997.7, 2014.3),
breaks = seq(1998, 2014, by=2), expand = c(0.02,0.02))+
scale_color_manual(values=my_colours) +
scale_fill_manual(values=my_colours)
Model4b_plot
#ggsave(file="Model4b_plot.png", plot=Model4b_plot, width=4, height=4)
#ggsave(file="Model4b_plot.pdf", plot=Model4b_plot, width=4, height=4)
# 1. Build lme4 model:
model5b <- lmer(
Coral_cover ~ Thermal_regime * Year +
(1|Region), data=ENSO.83_14, REML = FALSE)
summary(model5b)
## Linear mixed model fit by maximum likelihood . t-tests use
## Satterthwaite's method [lmerModLmerTest]
## Formula: Coral_cover ~ Thermal_regime * Year + (1 | Region)
## Data: ENSO.83_14
##
## AIC BIC logLik deviance df.resid
## 1443.8 1475.1 -711.9 1423.8 159
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.52432 -0.48425 -0.05873 0.47939 3.06957
##
## Random effects:
## Groups Name Variance Std.Dev.
## Region (Intercept) 94.08 9.699
## Residual 234.31 15.307
## Number of obs: 169, groups: Region, 13
##
## Fixed effects:
## Estimate Std. Error df
## (Intercept) -1737.0757 506.0616 162.8735
## Thermal_regimeTropical upwelling 1854.0315 708.3227 167.0809
## Thermal_regimeEquatorial upwelling 291.8845 1024.3367 165.7040
## Thermal_regimeSeasonal 2738.6990 918.4800 167.5534
## Year 0.8799 0.2530 162.6103
## Thermal_regimeTropical upwelling:Year -0.9207 0.3539 166.9356
## Thermal_regimeEquatorial upwelling:Year -0.1475 0.5117 165.5744
## Thermal_regimeSeasonal:Year -1.3650 0.4588 167.3530
## t value Pr(>|t|)
## (Intercept) -3.433 0.000758 ***
## Thermal_regimeTropical upwelling 2.617 0.009671 **
## Thermal_regimeEquatorial upwelling 0.285 0.776038
## Thermal_regimeSeasonal 2.982 0.003293 **
## Year 3.478 0.000648 ***
## Thermal_regimeTropical upwelling:Year -2.601 0.010116 *
## Thermal_regimeEquatorial upwelling:Year -0.288 0.773525
## Thermal_regimeSeasonal:Year -2.975 0.003365 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) Thr_Tu Thr_Eu Thrm_S Year T_Tu:Y T_Eu:Y
## Thrml_rgmTu -0.714
## Thrml_rgmEu -0.494 0.353
## Thrml_rgmSs -0.537 0.383 0.265
## Year -1.000 0.714 0.494 0.536
## Thrml_rTu:Y 0.715 -1.000 -0.353 -0.383 -0.715
## Thrml_rEu:Y 0.494 -0.353 -1.000 -0.265 -0.494 0.353
## Thrml_rgS:Y 0.537 -0.383 -0.265 -1.000 -0.536 0.383 0.265
anova(model5b) # Thermal regime and its interaciton with year are significant
## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Thermal_regime 2922.20 974.07 3 167.16 4.1572 0.007168 **
## Year 580.02 580.02 1 166.90 2.4755 0.117530
## Thermal_regime:Year 2898.06 966.02 3 167.02 4.1228 0.007496 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
ranova(model5b) # ramdom effects are significant
## ANOVA-like table for random-effects: Single term deletions
##
## Model:
## Coral_cover ~ Thermal_regime + Year + (1 | Region) + Thermal_regime:Year
## npar logLik AIC LRT Df Pr(>Chisq)
## <none> 10 -711.89 1443.8
## (1 | Region) 9 -718.96 1455.9 14.146 1 0.0001692 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
drop1(model5b)
## Single term deletions using Satterthwaite's method:
##
## Model:
## Coral_cover ~ Thermal_regime * Year + (1 | Region)
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Thermal_regime:Year 2898.1 966.02 3 167.02 4.1228 0.007496 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# 2. Predict values:
pred5 <- predict(model5b,re.form = NA)
#3. Bootstrap CI:
boot5 <- bootMer(model5b, predict, nsim = 1000, re.form = NULL) # Include ramndom effects
std.err <- apply(boot5$t, 2, sd)
CI.lo_5 <- pred5 - std.err*1.96
CI.hi_5 <- pred5 + std.err*1.96
#Plot
Model5b_plot<- ggplot(ENSO.83_14, aes(
x = Year, y = Coral_cover, colour = Thermal_regime)) +
geom_line(aes(y = pred5),size=2) +
geom_point(aes(fill=factor(Thermal_regime)),
shape = 21, colour = "black", size = 2,
stroke = 0.3, alpha=0.5) +
geom_ribbon(aes(ymin = CI.lo_5, ymax = CI.hi_5),size=2,
alpha = 0.03, linetype = 0) +
scale_y_continuous("Live coral cover (%)",
limits = c(-20, 90),
breaks = seq(0, 90, by=20), expand = c(0,0))+
scale_x_continuous("", limits = c(1983, 2014.3),
breaks = seq(1983, 2014, by=2), expand = c(0.02,0.02))+
scale_color_manual(values=my_colours) +
scale_fill_manual(values=my_colours)
Model5b_plot
#ggsave(file="Model4b_plot.png", plot=Model4b_plot, width=4, height=4)
#ggsave(file="Model4b_plot.pdf", plot=Model4b_plot, width=4, height=4)